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The  purpose  of  this  research  is  to  develope  and/flight 
test  concepts  for  the  estimation  of  angle  of  attack^ Ca)  and 
sideslip  (($)  using  an  inertial  reference  platform.  This 
development  was  further  broken  down  into  real-time,  inflight 
estimation  of  a  and  post-flight  estimation  of  a  and  /3._ 
Following  theoretical  development,  the  concept?-  were  tested 
with  NASA  F-15A  flight  data  and  examined  real-time  during  a 
NASA  Highly  Integrated  Digital  Engine -Control  CHIDEC)  flight 
test  using  the  F-15A  aircraft. 

Angle  of  attack  is  a  critical  parameter  in  the 
maneuverable,  high  performance  aircraft  of  today.  Yet  many 
errors  are  present  in  the  current  methods  of  obtaining  this 
angle.  An  accurate  method  of  d  and  (3  estimation  could 
eliminate  the  need  for  such  probes,  and  allow  these 
quantities  to  be  used  for  a  broad  range  of  applications.  An 
inflight  estimator  was  developed  for  computational  speed  and 
accuracy  using  inertial  navigation  system  linear 
accelerations  and  angular  gates.  A  second  system  based  on 
linear  recursive  modeling  was  developed  for  post-flight 
estimation  of  a  and  The  data  and  programs  specified  in 

this  research  are  applicable  only  to  thos«r  aircraft 

f\  '  ,  • 

mentioned,  but  the  methods  of  estimation  are  universal,  j  i  / 

I  would  like  to  thank  my  thesis  advisors.  Dr.  Robert 
Calico  and  Maj.  Daniel  Gleason  for  their  help  and  direction 
in  this  project.  In  addition,  I  am  indebted  to  NASA  for 
allowing  me  the  use  of  their  facilities  and  test  aircraft  for 
this  work.  The  superb  efforts  of  NASA  flight  control 
engineer  Heather  Lambert  turned  a  theory  into  reality. 

ii 


i 


tt 


n  for 


17 


j 


Thanks  also  go  to  NASA  Test  Pilot  Jim  Smolka  for  flying 
my  "robustness  maneuver".  Finally,  I  would  like  to  dedicate 
this  thesis  to  my  parents  who  have  been  a  tremendous 
inspiration  to  me  in  my  Air  Force  career. 


£ 

ill 


TABLE  OF  CONTENTS 


PAGE 


Preface  . ii 

List  of  Figures . vi 

List  of  Tables  . viii 

List,  of  Symbols  and  Abbreviations  . ix 

Abstract  . xiii 

I.  Introduction  . . 1 

Overview . 1 

Problem . 2 

Solution . 3 

Scope . S 

II.  Background  . 7 

Previous  Approaches . 7 

Proposed  Approach  . 9 

III.  Inflight  Angle  of  Attack  Estimator  Theory  ...11 

Basic  Lift  Equation . 11 

Moment  Summation  . 12 

Load  Factor  Determination . 18 

No-Wind  Angle  of  Attack . 22 

Three-Dimensional  Windage  Calculation  .  .  23 
Flow  Diagram . 24 

IV.  Post-Flight  Estimation  Theory . 26 

Basic  Equations  of  Motion . 26 

Estimator  Theory . . . 32 

Programs  DKF  and  DKFLAT . 37 

Filter  Tuning . 37 

V.  Estimator  Validation  . 40 

Validation  Objective . 40 

Inflight  Estimator  Validation . 40 

Computer  Test  Programs . 43 


•  % 


TABLE  OF  CONTENTS  < CONCLUDED) 


Page 


Post-Flight.  Validation  Test  Runs . 45 

Results  of  Validation  Run . 47 

VI.  Flight  Test . 54 

Test  Philosophy . 54 

Test  Aircraft . 55 

Instrumentation . 55 

Data  Reduction . V . 59 

Test  Methods  and  Conditions . 59 

Results  and  Analysis . 62 

VII.  Maneuvering  Flight  Demonstration . 73 

Purpose . 73 

Scope . . . 73 

Maneuver . 73 

Results  and  Analysis . 74 

General  Comments . ' . 78 

VIII.  Conclusions  and  Recommendations . 80 

General  . 80 

Objective  1  . 80 

Objective  2  . 81 

Objective  3  . 81 

Objective  4  . 82 

Objective  5  . 83 

Appendix  A:  Computer  Program  DKF . 83 

Appendix  3:  Computer  Program  DKFLAT . 100 

Appendix  C:  Computer  Program  LONG . 103 

Appendix  D:  Computer  Program  LAT . 112 

Appendix  E:  F-15A  Aircraft  Description . 120 

Appendix  F:  Key  INS/Air  Data  Outputs . 124 

Appendix  G:  F-15A  Stability  Derivative  Models  . 128 

Appendix  H:  ELEXI  Program  A0A . 138 

Appendix  I:  Flight  Test  Card . 144 

Bibliography . 147 

Vitae . 150 

v 


jyw 


LIST  OF  FIGURES 


Figure 

1. 

2. 

3. 

4. 

5. 

6. 

7. 

8. 

9. 

10. 

11. 

12. 

13. 

14. 

15. 

16. 

17. 

18. 

19. 

20. 

21. 

22. 

23. 

24. 
23. 
El. 


Page 


Basic  Longitudinal  Aerodynamic  Forces . 11 

Moments  Acting  on  an  Aircraft  in  Flight . 13 

Relationship  of  Load  Factor  and 

Normal  Acceleration  . 18 

Inflight  Estimator  Flow  Diagram . 25 

Body  Axis  System  . 26 

Angle  of  Attack  as  a  Velocity  Ratio  . 27 

Sideslip  Angle  as  a  Velocity  Ratio  . . 30 

Covariance  Analysis  for  Kalman  Filter 

Tuning . 38 

Comparison  of  LONG  and  D5CF  Derived 

Angle  of  Attack  . 44 

Angle  of  Attack  Recovery  using  INS 

Covariance  Model  . 46 


Angle  of  Attack  Recovery  C . 5  degree  Elevator)  . 48 
Angle  of  Attack  Recovery  <2  degree  Elevator)  ..49 


Sideslip  Recovery  C.5  degree  Rudder)  . 50 

Sideslip  Recovery  C2  degree  Rudder)  . 51 

Sideslip  Recovery  <2  degree  Aileron)  . 52 

NASA  F-15A  HIDEC  in  Flight . 56 

F-15A  HIDEC  YAPS  Noseboom . 57 

Test  Point  1  Level  Acceleration  Results . 63 

Test  Point  2  Level  Acceleration  Results . 65 

Maximum  Power  Level  Acceleration  Results . 67 

Abrupt  Pitch  Maneuver  Results . 69 

Test  Point  7  -  3  g  Vind  Up  Turn  Results . 70 

Test  Point  8  -  5  g  Vind  Up  Turn  Results . 71 

Robustness  Maneuver  #1  . 73 

Robustness  Maneuver  #2  . 77 

F-15A  General  Aircraft  Layout . 121 

vi 


> •  v.v;av 


LIST  OF  FIGURES  (CONCLUDED) 


Figure  Page 

FI.  F-15A  HIDEC  INS  . 125 

F2.  F-15A  HIDEC  Air  Data  Computer  System . 126 

Gl.  F-15A  Subsonic  Trimmed  Lii  t  Curve  <U)  . 130 

G2.  F-15A  Three  Dimensional  a  Estimator  Model  ....131 

G3.  F-15A  Zero  Lift  Angle  of  Attack  and 

Pitching  Moment  Characteristics . 132 

G4.  F-15A  Gross  Weight  vs.  Center  of  Gravity  (U)  .133 

G5.  F-15A  Moments  of  Inertia  . . 134 

G6 .  F-15A  Estimator  IXX  Model  . 135 

G7 .  F-15A  Estimator  IYY  Model  . 136 

G8 .  F-15A  Estimator  IZZ  Model  . 137 

II.  NASA  F-15A  HIDEC  Flight  Test  Card  - 

Robustness  Test . 146 


LIST  OF  TABLES 


Table  Page 

I.  USAF  Standard  INS  Measument  Accuracies . 31 

II.  Static  Inflight  Angle  of  Attack  Test  Points  ....11 

III.  Static  Inflight  Estimator  Validation  . 12 

IV.  Post-Flight  Validation  Test  Points . 15 

V.  Level  Acceleration  Evaluation  . 60 

VI.  Abrupt  Pitch  Evaluation  . 61 

VII.  Vind-Up  Turn  Evaluation  . 61 

VIII.  Robustness  Evaluation  . 71 

FI.  Key  Sensor  Outputs . 127 


viii 


1 


LIST  OF  SYMBOLS  AND  ABBREVIATIONS 


Symbol 


Title 


Units 


A 

Stability  matrix 

a 

Normal  acceleration  in 

BZ 

body  axis 

feet/sec2 

B 

Control  matrix 

-  - 

c 

Mean  aerodynamic  chord  of  wing 

feet 

G 

Ving-body  coefficient  of  lift 

-  — 

wWB 

G 

Zero  lift  pitching  moment 

m 

o 

coefficient 

6 

acceleration  of  gravity 

feet/sec2 

H 

Matrix  of  available  measurements 

-  - 

h 

Altitude 

feet 

I 

Identity  matrix 

—  ~ 

Ix 

Moment  of  inertia  about  x  axis 

slug-feet2 

iy 

Moment  of  inertia  about  y  axis 

slug-feet2 

Iz 

Moment  of  inertia  about  z  axis 

slug-feet2 

Ixz 

Cross-product  of  inertia 

slug-feet2 

ICG 

Kalman  gain  matrix 

—  — 

L 

Total  aircraft  lift 

pounds 

lt 

Tail  contribution  to  lift 

pounds 

I. 

Lift  of  wing-body  combination 

pounds 

WB 

M 

Mach  number 

m 

Moment  about  y  axis 

C pitching  moment) 

foot-pounds 

» 

Zero-lift  pitching  moment 

foot-pounds 

O 

N 

Moment  about  z  axis 

C yawing  moment) 

ix 

foot-pounds 

V 


4 


* 


A 


■1 


A 

I". 

h, 

•* 

•f 


ix 


Symbol 


Title 


Units 


n 

P 

PM 

PP 

Q 

q 

R 

r 


T 


u 


V 


v 

V 

X 

XM 

XP 


X 


VB 


Load  factor  g 

Roll  rate  radians/sec 

Error  covariance  prior  to 
measurement  -  - 

Error  covariance  after 
measurement  ~  - 

Matrix  of  model  uncertainty  -  - 

Pitch  rate 

Matrix  of  measurement 
uncertainty 

Yaw  rate 

V.i  ng  area 

Error  range  or  twice 
the  given  accuracy 

Body  axis  velocity  in 
x  direction 

True  airspeed 
Airmass  velocity 
Body  axis  velocity 
Inertial  velocity 
Wind  velocity 
State  vector 

State  estimate  prior  to 
measurement 

State  estimate  after 
measurement 

Distance  from  center  of  gravity 
to  wing-body  aerodynamic 
center  feet 


radians/sec 

radians/sec 

feet 


feet/sec 

feet/sec 

feet/sec 

feet/sec 

feet/sec 

feet/sec 


X  Distance  from  center  of  gravity 

T  to  horizontal  tail  aerodynamic 

center  feet 

2  Measurement  vector  -  — 

x 


Symbol 

Title 

Units 

a 

Angle  rf  attack 

degrees 

aa 

Initial  angle  of  attack  guess 

degrees 

f3 

Angle  of  sideslip 

degrees 

6 

Control  input 

_  - 

a2 

Variance 

-  - 

a 

Pitch  angle 

degrees 

Roll  angle 

degrees 

Heading  angle 

degrees 

SUBSCRIPTS  AND  SUPERSCRIPTS 

V 

Transpose 

-  _ 

Trim  value 

-  - 

- 

Perturbation  value 

-  _ 

. 

Time  derivative 

-  - 

> 

Primed  lateral  derivative 

_  - 

X 

x  direction 

-  - 

y 

y  lirection 

-  - 

z 

z  direction 

-  - 

N 

North  direction 

-  - 

C 

East  direction 

-  - 

W 

Vest  direction 

—  — 

xi 


1  ■'Oy 

.  >> 

Symbol 

Title 

Units 

» 

• 

ABBREVIATIONS 

a 

a 

> 

* 

a 

Aileron 

-  - 

AO  A 

Angle  of  attack 

degrees 

• 

• 

CADG 

Central  air  data  computer 

-  - 

i 

i 

i 

C.  L. 

Chord  line  of  nine 

-  - 

1 

* 

* 

E 

Elevator 

— 

a 

i 

HIDEC 

Highly  integrated  digital  engine 
control 

l 

INS 

Inertial  navigation  system 

—  — 

ics 

Inertial  reference  system 

-  - 

€l 

Lat 

Lateral 

-  - 

Lone 

Longitudinal 

— 

[ 

Rd 

Rudder 

-  - 

[ 

RllS 

Root  mean  squared  error 

-  - 

l 

Rw 

Relative  wind 

_  _ 

f* 

s 

T 

Thrust 

-  - 

w 

u 

• 

X 

u 

Partial  derivative 
example:  <5X 

sn 

r 

r 

r 

OD 

At  time  X 

seconds 

* 

<X  -  1) 

At  time  step  nrior  to  X 

seconds 

i 

r  V 

r  *>, 
f  «> 

/■ 

s 

' 

• 

*! 

xii 

ABSTRACT 


> 


fs 

I 


v>‘ 


‘7- 


Throughout,  aviation  history,  pilots  and  engineers  have 
had  to  rely  on  mechanical  angle  of  attack  (a)  and  sideslip 
Cf3)  probes  to  determine  an  aircraft's  position  relative  to 
the  airmass.  Recent  advances  in  the  stability,  accuracy  and 
reliability  of  inertial  navigation  and  reference  systems  now 
allow  angle  of  attack  and  sideslip  information  to  be 
calculated  from  internal  aircral  .  systems  and  a  central  air 
data  computer.  Conflicting  requirements  for  inflight  angle 
of  attack  information  and  post-f light  angle  of  attack  and 
sideslip  data  reduction  demand  two  seperate  methods. 

Inflight  algorithms  require  fast,  accurate  angle  of  attack, 
with  no  assumptions  on  vertical  wind.  Post-flight  usage, 
however,  demands  great  accuracy  with  no  assumptions  on  either 
sideslip  or  vertical  windage.  From  the  aircraft  equations  of 
motion,  angle  of  attack  and  sideslip  algorithms  will  be 
developed,  with  velocity  and  rate  inputs  of  the  type  expected 
from  an  aircraft  central  air  data  computer  and  inertial 
navigation  system.  A  computer  program  will  then  be  developed 
to  validate  these  equations.  A  Kalman  filter  algorithm  will 
also  be  designed  to  aid  in  estimating  data  output  from  these 
sources. 

Flight  test  will  consist  of  two  parts.  Initally, 
velocity,  rate,  acceleration,  and  aircraft  position  signals 
will  be  picked  off  a  standard  inertial  reference  system  on  a 
NASA  F-1SA  aircraft. These  signals  will  be  processed  using  the 
algorithms  developed,  and  estimated  angle  of  attack  will  be 
output.  These  outputs  will  be  compared  with  those 

xiii 
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predicted  by  the  inflight  and  post-flight  algorithm  for 
accuracy  and  usability.  A  demonstration  will  then  be 
conducted  on  the  NASA  F-15A  to  determine  usability  and 
accuracy  of  inertially  derived  angle  of  attack  information 
in  a  highly  maneuvering  environment. 

Conclusions  will  include  boch  the  usability  and  accuracy 
of  inertially  derived  angle  of  attack  and  sideslip. 
Applications  for  accurate  and  reliable  angle  of  attack  and 
sideslip  are  many.  Three-dimensional  windage  can  be  readily 
predicted  as  a  result  of  this  research,  along  wi\.h  airborne 
wlndshear  detection  and  stall  warning  systems.  Elimination 
of  mechanical  angle  of  attack  and  sideslip  probes,  with  their 
inherent  inaccuracy,  failure  rates,  and  time  lag,  will  also 
allow  for  angle  of  attack  information  to  be  used  as  reliable 
feedback  in  automatic  aircraft  control  systems. 


I 
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ANGLE  OF  ATTACK  AND  SIDESLIP  ESTIMATION 
USING  AN 

INERTIAL  REFERENCE  PLATFORM 


I.  Introduction 


vy. 


OVERVIEW  The  orientation  of  an  aircraft's  velocity  vector 
relative  to  the  airmass  surrounding  it  can  be  described  by 
two  angles,  the  angle  of  attack  (alpha  or  a)  and  the  angle  of 
sideslip  (beta  or  /3).  Since  all  aerodynamic  forces  and 
moments  depend  on  these  two  angles,  they  are  critical 
parameters  in  stability  and  control  flight  testing.  Indeed, 
angle  of  attack,  and  to  a  lesser  extent  sideslip,  are 
extremely  important  parameters  to  pilots  of  high  performace 
aircraft  in  an  operational  environment.  Angle  of  attack  is 
continuously  monitored  in  air-to-air  maneuvering  and  during 
approach  to  landing  In  addition,  the  two  angles  are 
utilized  for  stall  inhibitor  and  alpha  limiter  subsystems  in 
advanced  flight  control  systems.  But  angle  of  attack  and 
sideslip  are  extremely  difficult  to  measure  precisely. 
Throughout  the  brief  history  of  aviation,  both  pilots  and 
engineers  have  had  to  rely  on  mechanical  angle  of  attack  and 
sideslip  devices  to  determine  those  critical  angles  in  the 
three— dimensi-  tal  airmass.  If  angle  of  attack  and  sideslip 
could  be  accur-  -»I y  estimated  both  in  the  dynamic  inflight 
environment,  and  during  post-flight  data  analysis  sessions,  a 
great  deal  of  difficulty  in  performing  calibrations,  creating 
correction  curves,  and  calculating  performance  derivatives 
could  be  eliminated.  An  accurate,  real-time  knowledge  of 
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angle  of  attack  and  sideslip  could  aid  the  pilot  in  a  variety 
of  ways.  Just  as  accurate  post-flight  estimations  of  both 
critical  angles  could  aid  the  flight  test  engineer. 


PROBLEM  Currently,  angle  of  attack  and  sideslip  are 
measured  by  external  probes  mounted  on  the  forward  fuselage 
of  operational  aircraft.  These  rotating  mechanical  probes 
have  slots  which  then  align  the  probes  to  the  local  airflow, 
thus  providing  an  estimate  of  the  true  angle  of  attack  of  the 
ning-body  combination  in  the  airmass.  Flight  testing,  to  a 
great  degree,  relies  on  pitot  probe  mounted  angle  of  attack 
and  sideslip  vanes.  The  mounting  of  probes  far  ahead  on  the 
fuselage  eliminates  some  of  the  error  due  to  local  flow 
effects,  but  many  problems  still  exist  with  external  angle  of 
attack  (AO A)  and  sideslip  probes  rnd  vanes.  First,  being 
away  from  the  center  of  gravity  of  the  aircraft,  the  rotating 
probes  are  subject  to  pitch,  yaw,  and  roll  motion-induced 
errors  which  vary  with  Mach  number  and  maneuver,  increasing 
the  difficulty  of  calibration  dramatically!  1: 41 .  Second, 
being  mechanical  pick-off  devices,  these  probes  have  a 
definite  time  lag  that  hampers  their  use  in  flight  controls 
and  autopilot  applications! 2: 4431 .  This  time  lag  can  produce 
extremely  undesirable  characteristics  in  most  flight- 
controls.  Third,  several  airplanes,  most  notably  F-16  and 
F-lll  aircraft  have  been  lost  due  to  failure  of  the  angle  of 
attack  probe  inflight!31.  This  has  been  due  to  the 
requirement  for  stall  inhibitor  systems  which  attempt  to 
prevent  the  aircraft  from  exceeding  a  specified  angle  of 
attack.  The  flight  controls  will  respond  to  a  high  AOA  with 
a  pitch  down  maneuver  to  reduce  the  sensed  extreme  angle  of 
attack.  If  both  probes  have  failed  in  flight,  the  aircraft 
could  pitch  down  for  no  reason,  other  than  the  false  AOA 
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derived  from  the  failed  probe[ 4: 1-51 j .  For  such  a  vital 
instrument  on  a  high  performace  aircraft,  single  prcbe 
failure  could  result  in  aircraft  loss.  Finally,  the  very 
accuracy  of  the  mechanical  probe  in  determining  angle  of 
attack  and  sideslip  is  difficult  to  determine.  Thacker  found 
errors  in  angle  of  attack  over  the  subsonic  range  to  be  on 
the  order  of  1.3  to  2  times  the  actual  value  of  alpha,  and 
sideslip  error  of  1.73  to  2.4  times  actual  beta  values  for 
the  USAF/CALSPAN  NT-33A,  varying  with  Mach  numberCi:21.  The 
F-16  flight  control  system  currently  calculates  angle  of 
attack  with  twin  double-slotted  probes  to  an  accuracy  of 
approximately  0.4  to  0.3  degrees  over  the  entire  flight 
regimeC31.  These  probe  errors  are  applicable  to  all  probed 
aircraft  in  the  subsonic  range  of  flight.  The  realm  of 
supersonic  flight  makes  greater  demands  on  the  use  of  probes 
and  vanes  to  determine  angles  of  attack  and  sideslip. 
Increased  performance  dictates  a  greater  requirement  for 
knowledge  of  these  angles.  However,  at  extremely  high 
speeds,  penalties  from  drag  and  aerodynamic  heating  require  a 
“cleaner'*  method  for  angle  of  attack  and  sideslip 
determination  without  reference  to  external  devices  in  the 
airflow.  Indeed,  the  increased  reliance  on  angle  of  attack 
and  sideslip  information  at  these  critical  speeds  and  flight 
conditions  places  a  strong  requirement  on  sensor  redundancy 
for  critical  flight  control  systems.  That  requirement  could 
be  fulfilled  by  such  an  internal  AOA  and  sideslip  estimation 
system. 


SOLUTION  Angle  of  attack  and  sideslip  can  De  determined 
with  a  high  degree  of  accuracy  from  inertial  reference 
systems  and  central  air  data  computers.  The  current 
^  generation  of  ring  laser  gyro  inertial  navigation  systems 

<INS)  and  inertial  reference  systems  <IRS>  have  accuracies  in 
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VJv  pitch,  roll  and  heading  angles  on  the  order  of  3.6  x  10 

radians<.032  deg. >  and  angular  rate  measurment  accuracies  on 
the  order  of  7.3  x  iO'^radians/secondC . 043  degrees/sec) 
15:86-891.  The  newest  generation  ring  laser  gyro  INS  units 
are  even  better,  with  accuracies  of  2.0  x  10~®radiansC61 . 
Normal  acceleration  can  be  measured  to  an  accuracy  of  2 
feet/sec2[S] .  Using  the  three  INS  rates,  accelerations,  and 
central  air  data  computer  (CADO  inputs,  the  angle  of  attack 
can  be  estimated  with  a  moment  summation  and  lift  model. 

This  lift  can  then  be  compared  to  the  current  aircraft  weight 
model  and  measured  load  factor.  An  angle  of  attack 
estimation  can  then  be  generated  based  on  the  required  load 
factor.  For  post-flight  evaluation,  the  INS  rates  and  Euler 
angles  can  be  input  to  a  extended  Caiman  estimator  for  ground 
reconstruction  of  angle  of  attack  and  sideslip. 

The  speed  of  the  digital  INS,  combined  with  its  accuracy 
and  reliability  provide  the  following  advantages.  Primarily, 
if  angle  of  attack  can  indeed  be  estimated  with  a  high  degree 
of  certainty,  the  speed  of  the  INS/computer  calculations 
implies  that  the  information  can  be  used  in  conjunction  with 
advanced  flight  control  systems  as  a  feedback  quantity.  Being 
internal  to  the  aircraft  and  extremely  accurate,  the  INS  will 
eliminate  local  flow  and  Mach  effects  that  must  now  be 
corrected  in  raw  external  probe  data.  Finally,  the 
reliability  of  the  current  INS  systems,  and  even  higher 
reliability  of  the  ring  laser  INS's,  adds  sensor  redundancy 
to  probe  derived  AOA  systems.  As  the  mean  time  between 
failures  <MTBF)  of  INS's  increases,  it  could  become  the 
primary  source  of  angle  of  attack  and  sideslip  data. 

A  post-flight  derived  angle  of  attack  and  sideslip  can 
be  used  readily  to  calibrate  probes  and  aid  in  the  estimation 
of  stability  derivatives  without  lengthy  trim  shots.  Alpha 
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and  beta  derived  immediately  following  flight,  from  INS  and 
CADC  data  can  save  valuable  flight  time,  effort  and 
substantial  computer  resources. 
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SCOPE  The  purpose  of  this  thesis  is  to  develop  an  accurate 
method  for  estimating  angle  of  attack  and  sideslip  using 
signals  available  from  a  standard  INS.  The  estimation 
algorithm  is  broken  into  two  parts,  the  inflight  and  the 
post-flight  algorithms.  The  inflight  algorithm  must 
accurately  <to  .5  degrees)  measure  angle  of  attack  in  real 
time,  suitable  for  use  in  a  flight  control  application.  An 
update  rate  of  30  to  60  cycles  per  second  will  allow  flight 
control  usage  of  the  derived  information.  This  dictates 
quick,  simple  and  accurate  calculations  from  INS  and  CADC 
inputs  only.  A  by-product  of  alpha  esitmates  should  be  the 
accurate  calculation  of  three-dimensional  winds.  The 
post-flight  algorithm  must  accurately  C to  0.23  degrees) 
estimate  AOA  and  sideslip  using  INS  angles  and  rates.  These 
desired  accuracies  are  based  on  200  percent  improvement  on 
current  uncorrected  probe  accuracies  of  2  degrees  over  the 
entire  flight  regime.  This  information  would  be  input  to  a 
linear  recursive  estimator  for  reconstruction  of  the  required 
data.  There  is  no  requirement,  however,  for  real  time  data. 
In  this  way,  complicated  estimation  algorithms  can  be  applied 
to  determine  the  most  accurate  angle  of  attack  and  sideslip, 
along  with  windage,  for  flight  test  interests. 

The  algorithms  will  be  developed  with  these  specific 
purposes  in  mind.  A  complete  flight  test  program  will 
determine  the  usefulness  and  accuracy  of  the  angle  of  attack 
and  sideslip  estimates.  This  testing  will  first,  evaluate 
estimation  in  straight  and  level  flight,  and  gentle  turning 
maneuvers.  This  will  be  followed  by  examination  of  the 

-  3  - 


!SW 


%> 


5 


algorithms  in  a  highly  Maneuvering  environnent  at.  increased 
AOA  and  load  factors.  The  results  will  be  compared  to  a 
computer  simulation  of  the  angle  of  attack  and  sideslip 
required  for  the  specific  maneuver,  as  well  as  probe  and  vane 
measured  alpha  and  beta. 

The  objectives  of  this  project  are: 

1.  Develop  a  post-flight  algorithm  for 
determination  of  angle  of  attack  and  sideslip 
with  reference  to  INS  and  CADC  outputs  only, 
making  no  assumptions  on  vertical  wind  or 
sideslip.  The  estimations  should  be  accurate  to 
.25  degrees.  Then,  apply  a  Kalman  filter 
to  post  flight  data  for  accurate  AOA 
and  sideslip. 

2.  Validate  by  computer  simulation  the  angle  of 
attack  and  sideslip  recovered  by  the  Kalman 
filter. 

3.  Develop  an  inflight  algorithm  to  determine 
angle  of  attack  with  reference  to  INS  and  CADC 
outputs  only,  assuming  sideslip  is  zero.  The 
data  must  be  real  time  and  accurate  to  . 3 
degrees. 

4.  Demonstrate  algorithms  with  flight  test  data. 

3.  Determine  through  flight  test  where  the 
algorithms  are  no  longer  valid  estimators. 
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II.  BACKGROUND 


PREVIOUS  APPROACHES  Four  basic  approaches  have  been  used  to 
estimate  angle  or  attack.  Freeman's  method  C71  in  If""*' 
utilized  accelerometers  only,  with  no  inertial  hardware  The 
additional  aircraft  dynamics  were  then  estimated  by  control 
surface  deflection  pickoff  devices  and  basic  aircraft 
equations  of  motion.  An  angle  of  attack  processor  then  took 
inputs  of  flight  condition,  accelerations,  and  surface 
position  to  estimate  angle  of  attack  for  the  required 
maneuver.  This  method  required  extensive  modeling  of 

stability  derivatives.  Tail  aerodynamics  were  needed  to  ! 

i 

implement  the  algorithm.  Surface  position  indicators  were  j 

used  to  formulate  the  model  of  the  aerodynamic  response  of 

V 

\  *>  the  tail  surfaces.  Pitch,  roll,  each  and  (3-loading 

limitations  were  placed  on  the  estimator  to  simplify  this 
complex  modeling  into  a  managable  system  of  equations  for 
processing.  This  method  relied  heavily  on  aerodynamic 
modeling  using  stability  derivatives  and  surface  position  to 
evaluate  the  maneuver  being  performed  and  estimate  an  angle 
of  attack  required.  While  accurate  to  approximately  .  3 
degrees  angle  of  attack  over  the  flight  range  specified  in 
the  limitations  C7:383,  INS  information  could  eliminate  the 
need  for  such  extensive  aero  modeling  and  surface  position 
pickoff,  with  improved  accuracy  over  a  less  constrained 
flight  envelope. 

A  second  algorithm,  requiring  less  extensive  aerodynamic 
modeling,  is  provided  by  Petrov  and  Studnev,  et  al  (83.  This 
method  requires  precise  modeling  of  the  coefficients  of  lift  and 

At  drag  and  accelerometers  to  model  aerodynamic  performance  in 

'v’<y 
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gliding  flight  (8:31.  The  linearized  equations  are 
applicable  only  for  small  alpha  and  beta  (less  than  10 
degrees)  and  do  not  provide  the  angJes  throughout  the 
■aneuvering  envelope  [8:3]. 

Perhaps  the  most  detailed  work  has  been  accospl lathed  by 
Olhausen,  using  INS  outputs  for  YF-16  flight  test  £9].  The 
method  uses  INS  accelerometer  and  velocity  naasurments,  along 
with  Euler  angle  measuments.  Basic  aircraft  equations  of 
wiotion  were  solved  using  appropriate  order  Bunge— Eutta 
integration  techniques.  This  method  is  effective,  but  like 
most,  it  requires  the  assumption  that  sideslip  and  vertical 
windage  are  zero.  Used  primarily  for  flight  test 
applications,  the  algorithm  develops  problems  determining 
winds  in  steep  turns  where  sideslip  is  not  negligible. 

The  fourth  method,  and  the  one  used  by  Thacker  Ill,  is 

A 

state— space  estimation.  Thacker  only  used  two  states,  9  and 

A  A  A 

q  .  In  this  model ,  9  represents  the  pitch  angle  and  q  is  the 
pitch  rate.  These  quantities  are  perturbation  values  about 
some  nominal  flight  condition.  The  two  states  were  shown  by 
Thacker  to  successfully  determine  angle  of  attack  to 
approximately  .3  degrees.  Logically,  a  more  accurate  math 
model  of  the  aircraft  dynamics  should  result  in  more  accurate 
determinations  of  angle  of  attack.  In  addition,  a  model  of 

the  lateral  dynamics  of  the  same  aircraft  should  also  yield 

comparable  results  with  sideslip.  The  usual  model  of  an 
unaugmented  aircraft  consists  of  4  states  in  both  the 
decoupled  longitudinal  and  lateral  modes.  Although  varied, 

A 

this  mathematically  more  accurate  system  can  be  one  with  V, 

*  A  A  A  A 

a  ,  q  ,  and  9  ,  where  a  is  angle  of  attack,  V  is  true 

A  A 

airspeed,  0  is  pitch  angle,  and  q  is  pitch  rate.  Again, 
these  quantities  are  perturbation  values.  They  are,  in  other 
words,  the  changes  in  those  variables  from  some  steady  state 


values.  Accurate  inputs  can  be  obtained  for  q  and  6  from 
the  INS,  with  a  moderately  accurate  V  from  the  CADC  Coff  by 
the  angle  of  attack  rotation).  Alpha  would  then  be  estimated 
by  a  Kalman  estimator.  This  state— space  estimator  has 
several  disadvantages.  The  system  must  model  about  34 
stability  derivatives,  depending  on  the  accuracy  desired  in 
the  result.  These  derivatives  are  naturally  based  on  Mach, 
altitude,  fligho  condition,  and  require  mathematical  changes 
for  flcp  deflectxon.  This  extensive  math  modeling,  along 
with  the  computational  time  required  to  solve  the  system  of 
equations,  currently  precludes  this  method  for  real  time 
inflight  estimator  applications. 


PROPOSED  APPROACH  Angle  of  attack  and  sideslip  estimation 
has  two  inherently  opposing  requirements,  speed  and  accuracy. 
Speed  of  calculation  is  a  critical  requirement  for  flight 
^  1  control  usage  of  alphe,  as  well  as  for  pilot  Information  and 

three  dimensional  windage  derivation.  Calculations  must  be 
minimized,  with  a  Judicious  use  of  assumptions,  while 
retaining  accuracy  to  .3  degrees.  On  the  other  hand,  flight 
test  standards  require  an  accurate  knowledge  of  alpha,  beta, 
and  three-dimensional  winds.  This  information  is  not  time 
critical,  and  can  be  post-flight  processed  for  analysis  at  a 
later  time.  These  requirements  drive  the  need  for  two  types 
of  estimators,  a  rapid  estimator  method  for  inflight  use,  and 
a  very  accurate  estimator  method  for  flight  test  analysts 
after  the  fact.  I  propose  to  develop  these  two  methods  of 
angle  of  attack  and  sideslip  estimation,  and  flight  test 
their  validity. 

Freeman's  work  with  the  alpha  estimator  provides  a  good 
basis  for  the  inflight,  real-time  portion  of  the  alpha 
estimator  problem.  Incorporation  of  INS  rates,  angles  and 
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accelerations  can  dramatically  reduce  the  requirement  for 
extensive  stability  derivative  modeling,  and  increase  the 
accuracy  and  calculation  speed.  However,  the  flight  envelope 
imposed  by  Freeman  must  be  expanded  for  a  wider  range  of 
inflight  applications.  A  new  algorithm  based  on  total  lift 
and  aircraft  moments  will  be  - eveloped.  No  assumptions  will 
be  made  concerning  the  existence  of  a  vertical  wind,  and 
accurate  three-dimensional  windage  output  will  be  a  criterion 
for  acceptable  operation  of  the  inflight  alpha  estimator. 

With  calculation  speed  not  a  factor  for  the  post-flight 
estimator,  the  state-space  me  .hod  provides  a  basis  for  alpha, 
beta,  and  three-dimensional  windage  estimation,  with  no 
assumptions  made  on  any  of  these  three  quantities.  The  more 
extensive,  4  component  state-space  model  of  34  derivatives 
for  decoupled  longitudinal  and  lateral  response  modes  can  be 
easily  used.  This  would  also  require  application  of  a  linear 
recursive  Kalman  estimator  to  take  the  limited  INS  and  CADC 
inputs  to  estimate  the  alpha  and  beta  obtained  in  the 
maneuvers  performed.  3-d  windage  profiles  can  then  readily 
be  calculated  for  test  engineering  usage. 


III.  INFLIGHT  ANGLE  OF  ATTACK  ESTIMATION  THEORY 

BASIC  LIFT  EQUATIONS  Angle  of  attack  is  directly  related  to 
the  coefficient  of  lift  of  the  wing-body  combination.  If  the 
total  lift  on  the  wing-body  is  obtained,  angle  of  attack  may 
be  subsequently  computed.  Total  lift  is  the  sum  of  the  lift 
of  the  wing-body  combination  acting  at  the  corresponding 
wing-body  aerodynamic  center,  and  the  lift  of  the  tail, 
acting  at  the  aerodynamic  center  of  the  empennage. 


Figure  1.  Basic  Longitudinal  Aerodynamic  Forces 


Denoting  this  lift  of  the  wing-body  combination  as  Lw  , 
lift  of  the  tail  as  LT,  the  total  lift,  L,  is  given  by: 


L  -  L  ♦  L 

VB  T 


But  total  lift  can  also  be  defined  by: 


L  ■  nV 


where  n  is  load  factor  and  V  is  total  aircraft  weight. 
Then: 


nV  -  L  ♦  L 
VB  T 


C3> 


Ving-body  lift  can  be  expressed  in  a  standard  non-dimensional 


coefficient  G,  ,  where: 
lvb 


C4> 


VB 


<*pV2)S 


where  p  is  atmospheric  density  at  flight  altitude,  V  is  true 
airspeed,  and  S  is  wing  area. 

Rearranging: 


L  -  G  <ipV2DS  C3> 

VB  SfB 

so  that  in  terms  of  load  factor: 

nW  -  C  <|pV2)S  +  L  C6> 

MOMENT  SUMMATION  Lift  on  the  tail  must  be  modeled  next. 

This  can  be  accomplished  by  applying  an  analagous  equation, 
but  the  effects  of  elevator  deflection,  downwash,  and  surface 
position  pickoffs,  as  well  as  wake  effects  must  be 
considered.  These  effects  are  difficult  to  account  for  due 
to  the  inaccuracies  of  determining  exact  tail  deflection,  as 
well  as  the  errors  in  aerodynamic  modeling  of  the  flow  over 
the  empennage,  and  the  resulting  forces.  While  it  is 
possible  to  model  downwash  and  analytically  determine  tail 
lift,  the  equations  quickly  become  unmanagable  in  even  slight 
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maneuvering.  In  addition,  the  algorithm  would  then  vary 
extensively  with  tall  surface  control  design.  The  simpler 
and  more  direct  method  is  the  calculation  of  the  lift 
contribution  of  the  tail  surfaces  as  determined  through 
aircraft  moments. 

The  moment  action  on  the  total  aircraft  can  be  described 
adequately  by  the  zero-lift  pitching  moment,  wing-body  lift, 
and  the  moments  due  to  the  tail  lift.  Since  the  lift  is 
always  perpendicular  to  the  relative  wind,  the  resulting 
moments  and  arms  are  shown  in  Figure  2. 


Figure  2.  Moments  Acting  on  an  Aircraft  in  Flight 

where  mQ  is  the  zero-lift  pitching  moment  about  the 
quarter-chord  and  XT  is  the  distance  between  the  center  of 
gravity  of  the  aircraft  and  the  aerodynamic  center  of  the 
horizontal  tail.  X__  is  the  corresponding  distance  between 
the  center  of  gravity  of  the  aircraft  and  the  wing-body 
aerodynamic  center.  It  is  assumed  that  the  weight  acts 
through  the  center  of  gravity,  and  is  therefore  not  a 
contributing  factor  in  any  moment  equation.  This  is  a  valid 
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assumption  in  that  the  aircraft  in  flight  will  generally 
rotate  about  its  center  of  gravity.  Therefore,  any  moment 
summation  about  the  center  of  gravity  eliminates  weight  of 
the  aircraft  as  a  contributing  factor  in  the  moment 
equation.  A  normal  static  margin  will  place  the  center  of 
gravity  less  than  ten  percent  of  the  mean  aerodynamic  chord 
away  from  the  aerodynamic  center.  A  further  assumption  is 
that  the  lift  of  the  wing-body  combination  acts  at  a  single 
point  called  the  wing-body  aerodynamic  center.  These  two 
basic  assumptions  are  generally  held  to  be  true  for 
conventional  aircraft  throughout  a  large  portion  of  their 
flight  regimes.  In  addition,  centerline  thrust  is  assumed, 
resulting  in  no  moments  due  to  thrust,  as  well  as  drag. 

This  assumption  is  valid  for  NT-33 A  and  other  fighter-  type 
aircraft  investigations,  but  may  require  correction  in 
tankerxtransport  applications. 

The  term  mQ  as  used  here  refers  to  the  aero-lift  pitching 

moment  of  the  wing-body-tail  combination.  As  such,  the 

conditions  at  zero  lift  require  that  lift  of  the  wing-body 

offset  the  lift  of  the  tail.  Thus  the  terms  L_  and  L 

t  vs 

actually  refer  to  incremental  lift  from  the  zero-lift 
condition.  However,  the  zero- lift  values  of  wing— body  and  tail 
lift  are  negligible  when  compared  to  total  wing— body  and  tail 
lift  in  1-g  flight.  The  simplification  will  therefore  be 
made  that  L  and  L,_  are  total  lift  terms  and  not  incremental 
terms  from  the  zero-lift  condition.  So  the  summation  of  moments 
yields: 


In 


n  -  LX  ♦  L  X 

O  T  T  VB  VB 


<75 


where  is  the  wing-body  static  margin  as  defined  by  the 
expression: 

XVB  "  <h  -  hVB>5  <8> 


^v.V.  A  .VAW.V.V; J'.v.v.  A 


-  14  - 


*//,  But.  a  standard  alternate  formula  for  the  total  sum  of  the 

moments  Is  given  in  McRuer  £2:220]. 


Zm  -  qIY  +  pr<Ix  -  V  -  r  1^ 


P2I 


xz 


CO) 


where  q  is  pitch  acceleration,  p  is  roll  rate  and  r  is  yaw 
rate,  all  quantities  available  directly  from  a  USAF  standard 
INS.  I.,,  I„,  and  I_  are  respective  principle  moments  of 
inertia  to  the  X,  Y,  and  Z  axes  and  I___  is  th*?  applicable 
cross  product  of  inertia  calculated  in  the  body  axis. 

If  we  neglect  the  cross  products  and  rate-squared  terms 
as  negligible,  equation  (8)  becomes: 


Zm  ■  qly  +  pr<Ix  ~  <10> 

This  is  a  valid  assumption  due  to  the  relative  size  of  the 
cross-products.  The  NT-33A  data  £13:221  shows  a  cross 
product  of  480  slug-feet2  compared  to  a  difference  in  Ix 
and  Ix  of  20,000  slug-feet2  .  The  yaw  rate  squared  term 
will  always  be  less  than  1  radian/second,  and  it  is  assumed 
that  roll  rates  will  be  less  than  1  radian/second  also.  A 
typical  transport  aircraft,  the  Convair  88QM  £13:200], 
specifies  this  cross-product  as  0,  indicating  that  it  can  be 
ignored  for  the  purposes  of  this  investigation.  Now,  only 
the  moments  of  inertia  need  to  be  modeled  throughout  the 
flight  regime. 

So,  combining  equations  <73  and  <93  yields: 


qly  «■  pr<Ix  -  Iz3  -  mQ  -  LTXT  +  L^X^ 
Now  solving  for  LT  : 

L  -  j**  ~  Pr<1*  ~  ^  **0*  LvbXw 


<113 


<123 
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This  equation  for  the  lift  of  the  tail  provides  several 
advantages.  The  algorithm  does  not  require  extensive 
aerodynaaic  modeling  of  tail  effects.  The  theory  is  based 
entirely  on  the  effect  of  the  developed  tail  lift  on  the 
aircraft  pitch,  roll  and  yaw  motions.  In  that  light,  there 
is  no  downwash  calculation  error  or  surface  position 
indicator  requirement.  Indeed,  the  algorithm  doesn't  care  if 
the  horizontal  tail  surface  is  a  conventionally  flapped 
elevator,  a  full-flying  stabilator  with  variable  trim  tabs, 
of  even  a  differential  stabilator  with  roll  control  mixing. 

Equation  (11)  can  new  be  introduced  into  equation  (6) 
with  the  result: 


nV-C^  (|pV*)S+ 


vs 


-q!Y  -  Pr<lx-Iz>  ♦  Cm0(|pVa)Sc  ♦  L^X 


va 


(13) 


or 


X  nW-X  C  (ipV2)S 

T  T  lvb  2 


X„" 

•  .  _  ^*1 

i«- 

X, 

♦ 

T 

L  J 

(14) 


where: 


(|pVa)Sc 


Solving  for  C, 

S/B 


(13) 


estimator. 


yields  the  primary  equation  for  the 


v* 


nVXT  +  qly  +  pr<Ix-Iz)  -  Gm  <|pV2JSH 

**  o 


C16) 


WB 


C*pV2)SX„ 


1  +  J2L ' 

■  XT  - 


It,  can  then  be  noted  that  angle  of  attack  is  a  function  of 


WB 


a  -  FnCC,  ,  M,  h) 
lvb 


C17) 


with  M  being  Mach  and  h  being  altitude. 

This  functional  relationship  can  be  developed  for  each 

specific  aircraft  to  be  estimated.  Mach  number  and  altitude 

are  direct  CADC  outputs,  while  C,  can  be  derived  from 

uvb 

equation  <16). 

The  ratio  of  X^  to  XT  is  on  the  order  of  .02  for  most 
conventionally  stable  aircraft.  Thus  the  entire  correction 
factor  for  the  moment  due  to  the  lift  of  the  wing— body 
combination  is  approximately  1.02  to  1.03,  yielding  a  2  to  3 
percent  error  if  neglected  entirely.  This  term  will  be 
kept,  however,  to  increase  in-flight  angle  of  attack 
estimation  accuracy,  but  will  be  assumed  to  be  a  constant. 
This  neglects  the  center  of  gravity  shift  as  fuel  is  burned 
or  stores  are  released.  But,  this  shift  from  an  assumed 
medium  center  of  gravity  Jill  present  an  a  estimate  error  of 
less  than  1  percent.  Judging  by  normal  static  margin  shifts 
of  conventional  aircraft  through  a  normal  flight  mission. 

Also,  the  change  in  coefficient  of  lift  with  a  change 
in  a  is  assumed  to  be  instantaneous,  with  no  associated 
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dynamics.  Thus  the  angle  of  attack  required  for  a  specific 
coefficent  of  lift  to  exist,  as  determined  by  aircraft 
accelerations  and  momenta,  say  be  calculated  at  any  point  in 
time. 


LOAD  FACTOR  DETERMINATION  A  required  input  to  equation  (16) 


to  determine  C, 


is  load  factor  Cn).  The  load  factor  is 


defined  as  the  ratio  of  the  acceleration  of  the  aircraft 
normal  to  the  flight  path  within  the  plane  of  symmetry  to  the 
acceleration  of  gravity.  INS  linear  acceleration  would  then 
yield  only  acceleration  normal  to  the  aircraft  within  the 
plane  of  symmetry.  Vith  zero  angle  of  attack,  the 
relationship  of  a^  and  n  in  level  flight  is: 


■  -Cn  -  l)g 


(18) 


t.l.ero  aj  is  normal  acceleration  to  the  aircraft  and  g  is 
the  acceleration  of  gravity  at  the  earth's  surface,  which  is 
assumed  to  be  a  constant  for  the  flight  conditions  examined. 


Figure  3.  Relationship  of  Load  Factor  and  Normal 
Acceleration 
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This  equation  can  be  reasonably  corrected  for  large 
angle  or  attack  and  flight  path  angle  and  the  existence  of 
roll  by  rotating  the  acceleration  vector  to  a  position 
normal  to  the  flight  path.  Figure  3  depicts  the  difference 
in  orientation  of  the  load  factor  and  normal  acceleration  to 
the  flight  path. 

In  straight  and  level  flight,  normal  acceleration  is  0, 
but  load  factor  is  1.  However,  at  a  climb  angle  of  90 
degrees,  this  load  factor,  while  normal  acceleration  is 
0,  becomes  0  also.  The  same  is  true  in  a  bank  of  90  degrees 
while  maintaining  a  steady  course.  Normal  acceleration  is  0, 
and  load  factor  is  also  0.  This  change  from  normal 
acceleration  to  load  factor,  accounting  for  large  pitch  and 
roll  angles  can  be  found  by  rotating  the  body  normal 
acceleration  to  flight  path  normal  acceleration,  ayj,  in 
the  wind  axis  system. 

So,  using  the  transformation  matrix  Lwb  from  body  to 
wind  frames: 


cosocos/3 

sin/3 

sinacos/3 
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However,  it  is  immediately  evident  that  prior  knowledge  of 

v/ 

angle  of  attack  and  sideslip  is  required  for  this  rotation. 
As  this  is  only  used  for  a  small  correction  factor 
to  obtain  an  accurate  load  factor,  n,  some  simple 
approximations  may  be  used  for  a  and  (i.  For  this 
case  only,  sideslip  is  assumed  0.  A  no— wind  angle  of  attack 
can  be  quickly  computed  using  direct  INS  data.  Vhile  not 
accurate  for  estimator  purposes,  this  approximate  a  can 
serve  very  well  to  correct  body  axis  accelerations  to 
approximate  fligh  path  accelerations  for  use  in  equation 
(19).  The  equations  for  quickly  obtaining  this  a are 
presented  in  the  next  section  on  No-wind  Angle  of  Attack. 
Using  this  aQ  and  assuming  zezo  sideslip  for  the 
transformation  matrix  purposes  only,  becomes: 


vb 


coeaa  0 

0  1 


sinaa 

0 


-sinaa  0 


cosaa 


(21) 


Then,  using  the  transformation  matrix,  normal  acceleration  in 
the  wind  axis  system  is: 

a„  — (sina  )a_  +  (cosa  )a_  (22) 

Vx  O  Bx  a  Bx 

But  normal  acceleration  in  the  wind  axis  say  also  be  written 

as: 
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aWa-  (component  of  g^Ws  ■  -ng 


(23) 
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So  to  find  the  component  of  gravity  in  the  z  direction  in  the 
wind  axis,  another  transformation  matrix  is  used  to  rotate  g 
from  the  vertical  reference  frame  to  the  wind  axis.  It 
should  be  noted  that  wind  Euler  angles  are  used  in  this 
transformation.  However,  for  the  purposes  of  this  estimate 
of  n,  body  axis  Euler  angles  from  the  INS  will  be  used. 
is: 
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4 
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vv 


cosdcosy/ 


cosdsiny 


-sind 


sin^sindcosy/  sin^sindsiny*  sin^cosd 
— cos^sinv'  +ccS(pcos^ 


cos^sindcosv/  cos^sindsinv  cos<^cosd 
+sin<^siny/  -sin^cosv 


(24) 


Heading  angle  is  not  a  factor  in  this  determination,  so  it  is 
assumed  to  be  0.  Again,  a  is  assumed  to  be  small  (below  20 
degrees)  and  /3  is  assumed  to  be  0,  so  d*dw  and  LwV 

becomes: 
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wv- 


cosd  0 

sin^sind  cos^ 
cos^mind  — sin^ 


-sind 

sindcosd 

cos^cosd 


(23) 


Then  the  component  of  g  in  the  z  direction  in  wind  axis  is: 
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> 
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<j 
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i 
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gVg  ■  (cos^cosd)g 


Substituting  equations  (22)  and  (26)  into  equation  (23) 


(26) 
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********** 


yields  the  following  equation  for  load  factor  determination. 


(cosaa)aBx“<sina0)aBx  ■  -ng  *  gccsOcosd 


<271 


where  aa  is  the  initial  angle  of  attack  guess,  9  and  <p  are 
pitch  and  roll  angles. 


Solving  for  n  in  terms  of  body  axis  accelerations  and 


Euler  angles  yields: 


-<cosa  )a„  +  <sinerlam 

a  Bx  a  Bx 


+  (coo6)(cos^) 


<281 


This  load  factor  can,  derived  from  INS  accelerations  and 
Euler  angles,  become  an  input  to  the  primary  estimator 
equation  <161.  a0,  or  the  approximate  angle  of  attack,  is 

the  only  value  that  needs  to  be  calculated  for  use  In  the 
load  factor  equation  <281.  A  no-wind  assumption  will  allow 
quick  computation  of  this  rough  guess. 


NO-VIND  ANGLE  OF  ATTACK  A  no-wind  angle  of  attack  can  be 
immediately  obtained  from  the  INS.  All  inertial  velocities 
are  actually  groundspeeds  over  the  locally  flat  earth.  If 
we  assume  for  this  guess  only  that  wind  is  negligible  as 
compared  to  groundspeed,  the  INS  groundspeeds  will  also  be 
true  airspeeds,  and  an  initial  aQ  can  be  derived  from 
inertial  velocities  in  the  X,  Y,  and  Z  directions. 


vBj|  ■  <cosO  cosf*lv1N  +  <cos6  sinplvZB  -  <sin61vzz 


<291 


•  <-cosd  sinyr  ♦  Bind  sint>  cos*lv] 
< co*d  cosy  ♦  sind  sinO  sinylvI( 
<sind  cosGlv 


<301 
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v  "  (sirup  sirup  +  cos<p  sind  cos^)»XN 
C-sin^  cosy/  +  cos<p  sind  siny/)vIls 
(COS(p  COS03V 


These  body  velocities  can  then  be  evaluated  to  find  an 
approximate  angle  or  attack  using  the  body  axis  relationship: 


aQ  m  arctan 


This  no-wind  guess  of  angle  of  attack  is  then  used  only 
to  make  the  small  rotational  correction  on  normal 
acceleration  to  the  load  f actoz . 


THREE  DIMENSIONAL  VIN9ACE  CALCULATION  With  angle  of 


attack  recovered  from  the  functional  equation  (17),  airmasm 
velocities  of  the  aircraft  can  be  calculated  in  north,  east 
and  down  axes. 


(cosy/  cos0  cosa  cos/3  +  cosy/  sind  wirup  sin/3 
cosy/  sind  cos^  sina  '•os /3  -  siny/  cos<p  sin/3 

siny/  simp  sina  cos/3)  V 


v4_  ■  (siny/  coed  cosa  cos/3  +  siny*  sind  sin^  sin/3 

AS 

siny/  sind  cosp  sina  cos/3  +  cosy/  cosp  sin/3 
cosy/  sin^  sina  cos/3)  V 


v^  ■  (-sind  cosa  cos/3  +  cosd  sin<p  sin/3  + 


cosd  cosy  sina  cos/3)V 


where  v  ,  v  ,  and  v^  are  air 


velocities  of  the 


aircraft  in  the  direction  specified  and  V  is  the  true 
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airspeed  from  the  CADC. 


Assuming  sideslip  is  zero,  equations  (33),  (34),  and 
(33)  become: 

v  -  (cosy  cosd  coso  ♦  l3S y  sind  cosy  sina  +  (36) 

AN 

siny  siny  sina)V 

v  _  ■  (siny  cos©  coso  siny  sind  cosy  sina  -  (37) 

cosy  siny  sino)V 

v  _  ■  ( -oin©  coso  ♦  cos©  cosy  sina)V  (38) 

The  INS  groundspeeds  subtracted  from  these  airmass 
velocities  then  give  the  wind  component  in  each  direction. 
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FLOW  DIAGkjVM  The  overall  components  of  the  inflight 
e*  Imator  would  act  together  to  first  compute  the  no— wind 
guess  of  angle  of  atta  k.  At  this  point,  all  the  information 
necessary  for  est imation  of  act ual  angle  of  attack  would  be 
available  from  the  INS,  the  CADC,  and  the  model  of  the 
current  aircraft  weight  and  configuration.  A  CADC  output 
would  be  used  tc  determine  the  values  for  the  3  stability 
derivatives  that  must  be  modeled  for  varying  flight 
conditions.  Ti^e  angle  of  attack  would  then  be  computed  by 
the  estimator.  This  angle  of  attack  and  the  INS  groundspeeds 
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could  then  be  used  to  find  a  three— dimensional  wind,  given 
the  current  true  airspeed.  The  flow  diagram  in  figure  4 
depicts  this  action. 


Figure  4.  Inflight  Estimator  Flow  Diagram 
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IV.  POST  FLIGHT  ESTIMATOR  THEORY 


BASIC  EQUATIONS  OF  MOTION  The  body  axis  syste«  provides  a 
convenient  system  for  measurement  of  perturbed  ancle  of 
attack  and  sideslip  from  a  trimmed  condition.  The 
conventions  of  this  body  system  are  presented  in  figure  3. 


Figure  5.  Body  Axis  System 

In  the  body  axis  system,  the  angle  of  attack  can  be  related 
by  the  ratio  of  vertical  velocity,  v,  to  longitudinal 
velocity,  u. 


'  V 


Figure  6.  Angle  of  Attack  as  a  Velocity  Ratio 


Assuming  a  perturbation  away  from  the  trimmed  condition, 
this  ratio  takes  the  form: 


a  »  arctan 


(421 


In  order  to  determine  these  velocities  during  maneuvers, 
perturbation  equations  of  motion  can  be  developed  to 
simulated  the  interaction  of  pitch  rates  and  angles  on  the 
velocities,  U  and  V.  Development  of  dimensional  perturation 
equations  of  motion  can  be  found  in  McRuer(21.  Four 
linearized  equations  of  motion  can  be  derived  for  use  in  an 
angle  of  attack  estimator  algorithm.  The  first  of  these  is 
the  longitudinal  velocity  perturbation  equation. 

a  *  *  A  A  a 

u«*x  u  +  xw  +  x.w  +  Cx  -w  lq  -g  cos6_(0>  +x*_<5E  +  x ,  <5T  C431 
u  v  v  q  o  ^  •  o  ok  Or 


where  the  dimensional  derivatives  are  described  in  the 
glossary  of  terms.  Assuming  x^  to  be  small,  the  u  equation 
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i rJ  A  A  •*.  A  A  **,  A  •*,  /,  A  a  »r.  /,  r  t  a  a  *-  a  .•  •-  /.  a  <* 


k.v. 


can  be  written  as: 


u«xu  +  xw  +  Cx  -w)q  -g  cos©  C0)  *  x-<5E  ♦  x*_OT  C44) 

u  w  q  O  ^  O  OS  OT 


The  nornal  velocity  equation  is  fcrned  in  a  like  manner. 


w-z.w-z  u  ♦  z  w  *  Cz  +u  )q  -  r  sin©  C0)  ♦  z,_6E  ♦  z,r_OT  C45) 

wu  w  q  O  ^  “  O  OS  OT 


This  can  be  subsequently  reduced  to: 


z  u  z  w  Cz  +u,)q  g  sinC©„)0 

4.  v  4  q  °  _  o 

Ci-z.)  Cl-z.)  Ci-z. )  Cl-z. ) 

WWW  V 


*0s6E  ^  zteSE 


Cl-z.)  Ci-z. ) 

W  V 


C46) 


And,  in  the  saae  manner,  pitch  acceleration  is: 


q«au  +  a*+aq  +  a.w*  a.  6€  ♦  a  -  OT 

^  u  w  ^  ^  w  Ott  OT 


C47) 


where  w  is  defined  by  equation  C43).  Substituting  equation 
C43)  into  equation  C47)  gives: 


q  “  » 


a.  z  "I-  m.z  -  a.  Cz  -Mi  >1. 

*  v  u  u  ♦  av+  w  v  w  ♦  a  ♦  *  ««  0  q  - 

Cl-z.  )  V  Ci-z.)  *  Ci-z.) 

w  w  w 

J  It 


a.  g  sin©  1*  a.  z,_  a.  z  . 

v  0  ©  +  a  +  v  **  OE  ♦  a  ♦  v  ZL 

Os  Or 

Cl-z.)  Cl-z.)  Cl-z.) 

V  4  b  W  «■  W 


OT  C48) 
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A  fourth  equation  simply  supplies  the  identity: 


a  ■  q 


<450 


Equations  (44),  (46),  (48),  and  (49)  can  now  be  compiled 
into  standard  state-spac.2  matrix  notation  of: 


X  -  A-X  ♦  B  U 


(SO) 


Filling  in  the  stability  derivatives  gives  the  form  of  the  A 
and  B  matrices. 
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(31) 
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'V  In  an  analagous  nanner  to  the  development  of 

longitudinal  perturbation  equations ,  lateral  perturbation 
equations  can  be  defined  to  estimate  sideslip.  Beta  can  be 
defined  also  by  a  ratio  of  velocities,  in  this  case,  lateral 
velocity  to  true  airspeed. 


Figure  7.  Sideslip  as  a  Velocity  Ratio 
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<36> 


V  -  Yvv  +  wQp  +  g  cosC0n)^  -  u„r  ♦  Y^6A  +Y^«5Rd 


o  :5a  ^*<5Rd‘ 
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Thus  lateral  perturbation  equations  can  be  resolved  to 
the  following  four  equations: 


v  -  Y^v  -  Yvv  +  <Y  +wQ)p  ♦  <Yr-u0)r  ♦  g  cosO^* 


Y6a6A  *  Y<5*  4*** 


CS3> 


P  *  L„v  +  L„r  +  L;>v  +  L^dA  +  L^^Rd 


'ARd 


<54> 


<P  *  P 


(55) 


Assuming  that  Y^-O  ,  Y^-0  ,  Yr*0  ,  L^*0  ,  and  that 

P 

N^-0,  the  above  set  of  lateral  equations  reduce  to: 


P  *  L~v  ♦  KP  *  L.r  ♦  L6a6A  ♦  L^6Rd 


'<5Rdl 


<385 


4> 


(39) 


N  v 
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N  p  +  N  r 

P  r 


N<5a6A  +  ^d6^ 


(60) 


v5v 

vV 


v> 


This  in  matrix  notation  described  by  equation  (50)  is: 


Y  w  acoifl, 
v  O  O  O 


y<5a  Y«5Rd 


L-*  1-'  O  L.' 

v  p  r 


U*A  L6Rd 


o  o 


(61) 


N'  N'  O  N' 

v  p  r 


n6a  N6«d 


These  stability  derivatives  for  both  longitudinal  and  lateral 
A  and  B  matrices  can  be  evaluated  for  the  flight  condition  of 
interest.  The  numeric  A  and  B  matrices  can  then  be 
discretized  for  use  in  a  discrete  Kalman  estimator  program. 

ESTIMATOR  THEORY  Both  state-space  models  (31)  and  (61)  will 
approximate  their  respective  longitudinal  and  lateral  systems 
in  perturbed  motions.  Implied  in  this  is  the  assumption  that 
the  systems  are  linear  and  excited  by  small  perturbations 
about  a  trim  condition.  Equation  (30)  describes  a  state 
model  completely,  provided  all  states  are  available  for 
measurement.  This  single  equation  system  is: 


X  -  AX  BU 


(30) 


However,  all  the  states  of  any  given  system  may  not  be 
available  for  measurement .  In  this  case,  the  model  consists 
of  2  equations,  with  the  second  describing  the  available 


measurements.  This  system  is  formed  as: 


X  -  A-X  ♦  B  U 

A  A 

Z  -  H  ■  X 


(62) 

(63) 
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where  Z  is  the  measurement  vector  and  H  is  the  matrix  or  the 
available  measurement  components.  In  both  the  longitudinal 
and  lateral  cases  studied  here,  not  all  the  states  are 
available  for  measurement.  In  the  longitudinal  case,  only  Q 
and  q  are  available  from  the  INS.  In  the  lateral  case, 
the  only  measurements  are  r,  ,  and  p.  This  then  dictates 
the  H  matrix  of  available  measurement  components  for  use  in 
equation  (63). 


(64) 


(63) 


A  second  alteration  to  the  basic  system  equation  (30)  is 
the  presence  of  noise  in  both  the  modeling  process  and  the 
measuring  process.  The  modeling  cannot  take  into  account 
every  condition  and  outside  action  that  may  possibly  affect 
the  states.  Also,  the  model  itself  may  not  completely  depict 
the  exact  action  of  * ve  states  in  response  to  a  specified 
input.  This  inaccuracy  inherent  in  the  model  itself  can  be 


accounted  for  by  adding  a  proc 


noise.  The  errors  are 


random,  an<f  centered  about  a  zero  mean  of  the  actual  values. 
Thus  these  errors  can  be  represented  by  a  gaussian  noise 
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system  which  will  be  called  V  .  This  noise  will  be  defined 

by  a  variance,  o2  ,  with  a  mean  of  0.0.  In  additor.,  the 
available  measurements  are  also  corrupted  by  noise.  This 
noise  is  uncorrelated,  gaussian,  with  zero  mean.  This  noise 

V  can  be  described  by  a  variance,  o2  ,  analagous  to  the 
process  noise  describe  previously.  The  measurment  noise, 
affecting  the  measuring  device  itself,  affects  the  accuracy 
and  consistency  of  the  measurement.  Table  1  depicts  the 
accuracies  of  a  USAF  standard  INS.  The  jitter  value  is  added 
to  the  accuracy  to  produce  an  error  range,  T.  It  then 

follows  that  the  variance  of  the  Gausssian  approximation,  o2, 

is  j2<2T>2  for  use  in  the  measurement  uncertainty  matrix, 

R. 


Table  I. 

USAF  Standard  INS  Measurement  Accuracies 


SIGNAL 

UNIT 

ACCURACY 

JITTER 

e 

rad 

. 00028 

.0001 

* 

rad 

. 00028 

.0001 

p 

rad/sec 

. 00073 

. 00033 

q 

rad/sec 

. 00073 

. 00033 

r 

rad/sec 

.00073 

. 00033 

Assuming  that  the  matrices  Q  and  R  are  stationary,  the 


resulting  R  matrix  is: 


.ONO 


4.81  x  10 


4.03  x  10"7 


C66) 


LAT 


4.02  x  10-7 


4.81  x  10' 


0 

0 

4.03  x  10 


-7 


C67> 


A  similar  Q  matrix  contains  the  uncertainty  values  for 
the  model  itself.  These  covariances  will  be  the  subject  of 
later  filter  tuning  requirements. 

The  Kalman  estimator  makes  initial  estimates  of  the 

*  # 

states  prior  to  the  measurement,  and  then  subsequent  estimates 
after  the  measurement.  These  estimates  are  provided  by  the 
initial  conditions  placed  on  the  states,  followed  by 
subsequent  values  of  XP,  the  Kalman  predicted  state  vector. 

The  full  derivation  of  the  Kalman  equations  can  be  found 
in  GelbC121.  The  prediction  equation  takes  the  form: 


1 

*1 


1 


XMCK3  -  A  - XPCK—i}  +  B-IKK-O  C68) 

where  XM  is  the  state  estimate  prior  to  the  measument. 

In  a  likewise  manner,  the  error  covariance  prior  to  the 
measurwient  can  be  calculated,  using  an  initial  best  estimate 
of  covariance,  by: . 


I 

7, 


•  *>  v* 


where  PM  is  the  error  covariance  prior  to  the  Measurement  and 
PP  is  the  error  covariance  after  the  measurement,  at  time 
CK  -  1). 

With  these  estimates,  the  Kalman  gain,  KG,  to  be  applied 
to  the  post-measurement  error  covariances  and  state 
estimation  can  be  found.  This  Kalman  gain  matrix  is 
calculated  from  the  expression: 


KG(K)  -  PM<K3HT[  H.PM<K)  Ht  +  R  ] 


(70) 


The  predicted  error  covariance,  given  the  next 
measurement,  is  calculated  as: 


PPCK3  ■  [  I  -  KGC K) • H  ]?MCK-i> 


<71  > 


and  the  corresponding  state  estimate  is  calculated  as: 


XPCK3  -  XMCK3 


♦  KG  CIO  [ 


ZCK)  -  H-XMOO 


Now  given  a  full  state  estimate  in  both  longitudinal  and 
lateral  modes,  a  and  (i  are  estimated  from  the  perturbation 
and  true  velocities  by: 


a<K)  *  arc TAN 


w_  «■  w<K) 


uQ  ♦  u<K3 


C73> 


(KID"  arcSIN 


v0  +  v<K3 


<743 
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rhere  w„  ,u  .and  v„  are  trim  values  and  V  is  the  CADC 
o  ’  o  '  o 

calculated  true  airspeed. 

PROGRAMS  DKF  AND  DKFLAT  A  discrete  Kalman  filter  program  by 
Gleason  [131  was  used  to  perform  the  above  calculations. 

This  program  was  modified  for  angle  of  attack  and  sideslip 
determination  into  two  separate,  parallel  programs,  DKF  and 
DKFLAT.  A  sample  listing  of  DKF  is  given  in  appendix  A. 
DKFLAT  is  identical  to  DKF  with  the  exception  of  the  sideslip 
subroutine  BETA,  which  is  listed  in  appendix  B.  The  programs 
take  true  angle  of  attack  and  true  airspeed,  and  the 
discretized  A  matrix  for  lateral  and  longitudinal  modes. 

Step  elevator  or  thrust  changes  are  input  for  the 
longitudinal  responses  in  DKF,  resulting  in  a  trace  of 
computed  angle  of  attack,  and  estimated  angle  of  attack  using 
only  the  two  measurements  of  Q  and  q  .  Step  rudder  and 
aileron  inputs  are  inserted  to  DKFLAT,  resulting  in  a  trace 
of  computed  (3  and  the  estimated  (3  using  the  three 
measurements  of  r,  and  p. 

FILTER  TUNING  The  0  matrix  can  now  be  modified  to  provide 
the  proper  degree  of  uncertainty  to  the  model.  A  good 
starting  point  for  Q  seems  to  be  the  same  uncerta  ty  values 
as  in  R.  At  this  point,  a  covariance  analysis  can  be 
performed  to  match  computed  root  mean  squared  <RMS>  error  to 
the  true  RMS  of  the  system.  This  will  provide  the  optimal 
estimator  for  use  with  actual  flight  test  data. 

The  true  RMS  is  the  difference  between  the  filter 
estimate  and  the  true  value  of  the  system. 

The  difference  between  the  two  quantities  is  due  to  the 
weighting  given  to  the  measurements  as  opposed  to  the  model. 


Judged  by  model  formulation  and  sensor  capability.  These 
weightings  can  be  varied  through  choice  of  R  and  Q  matrices, 
based  on  known  measurement  errors  and  modeling 


Figure  8.  Covariance  Analysis  for  Caiman  Filter  Tuning 

inadequacies.  R  is  fixed  by  the  physical  constraints  of  the 
given  INS  accuracies  and  Jitters.  Q  can  be  adjusted,  though, 
until  approximate  RMS  equality  is  obtained.  Figure  8a  shows 
the  result  of  low  measurement  weighting,  while  8b  depicts  too 


much  weighting  on  the  actual  measurement.  Figure  8c  depicts, 
therefore,  a  tuned  Kalman  estimator,  with  correct  weightings 
applied  to  the  measurement  versus  the  model.  The  true  RMS  is 
approximately  equal  to  the  computed  RMS  in  this  tuned  case. 

At  this  point,  the  Kalman  estimator  is  ready  for  test  runs  to 


verify  its  operation. 


V.  ESTIMATOR  VALIDATION 


VALIDATION  OBJECTIVE  The  basic  objective  or  the  estimator 
validation  phase  is  to  ensure  that  both  estimators  have  been 
formulated  correctly,  with  a  minimum  error  under  known, 
static  conditions.  The  aerodynamic  data  available  will  allow 
a  study  of  the  effects  of  assumptions  and  neglected  higher 
order  terms  on  the  overall  accuracy  of  both  estimator 
systems,  in  the  absence  of  noise  and  jitter  associated  with 
actual  systems.  A  comparison  of  the  estimated  angle  of 
attack  and  sideslip  to  actual  alpha  and  beta  is  the  overall 
goal  of  the  validation.  Once  both  algorithms  have  been 
verified,  a  simulation  will  be  accomplished  to  determine  the 
estimator  characteristics  in  a  dynamic  environment,  close  to 
actual  flight  conditions.  The  simulator-derived  angle  of 
attack  and  sideslip  can  also  be  compared  to  estimated  values 
to  study  the  impact  of  neglected  dynamic  effects. 


INFLIGHT  ESTIMATOR  VALIDATION  Six  inflight  points  were  used 
to  verify  the  ability  of  the  inflight  estimator  to  recover  a 
and  f)  under  stable  flight  conditions.  NT-33A  data  from  the 
National  Aeronautics  and  Space  Administration  CNASA2 
contained  in  a  handling  qualities  report  CIS]  was  the 
reference  data  for  this  phase  of  the  validation. 

The  NT-33A  aircraft  is  a  programmable,  variab.'  e 
stability  aircraft  modified  from  a  basic  T-33  jet  trainer. 
This  two-place  aircraft  is  capable  of  a  wide  range  of  flight 
conditions,  and  a  great  deal  of  aerodynamic  data  has  been 
recorded  for  use  in  studies  of  this  type. 

In  a  clean  aircraft  configuration  with  a  nominal  flight 

,s'\ 

V>  control  system,  6  data  points  were  established,  as  shown  in 


Table  II,  to  relate  a,  altitude,  Mach,  and  angle  of  attack. 


Table  II.  Static  Inflight  Angle  of  Attack  Test  Points 

NT-33A  Aircraft 


DATA  POINT 

3 

4 

5 

6 

mm 

8 

247 

819 

62.  7 

222 

440 

129 

cu 

WB 

.2539 

.0891 

.9477 

.2805 

.  1504 

.4698 

MACH 

.4 

.7 

.3 

to 

to 

.75 

.65 

ALTITUDE(FT) 

0 

0 

20k 

20k 

20k 

40k 

a( RADS) 
TRUE 

.016 

-.016 

.164 

.014 

-.005 

SSB 

A  multiple  linear  regression  analysis  for  2  independent 

variables  was  performed  to  show  angle  of  attack  as  a  function 

of  C  and  Mach.  The  altitude  variance  was  not  sufficient  to 
>.  lwb 

provide  a  good  correlation  for  its  incorporation  into  the 
regression  formula.  This  analysis  was  accomplished  by  a 
least  squares  method  using  the  HP-41 C  hand  computer.  The 
regression  for  the  NT-33A  in  clean  configuration  provided  the 
following  formula: 


a(RAD)  -  -.02  +  . 2CC  )  -  . 03CM)  (73) 

uvs 


For  each  of  the  6  flight  points,  the  C  required  for 

UVB 

flight  was  found  through  the  primary  estimator  equation  (16). 
In  straight  and  level  flight  for  each  of  the  points,  the  load 
factor  is  1,  by  definition,  with  no  pitch,  roll  or  yaw  rates 
or  accelerations.  Equation  (16)  then  reduces  to: 
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C  _( |pV2>Sc 

mO  2 


(76) 


-v>‘ 


vx_ 


WB 


(|pV2)SX„ 

7.  T 


1  ♦  WB 


The  Ct  required  for  straight  and  level  flight  is  calculated 
lvb 

by  this  equation.  This  is  then  input,  to  the  regression 
formula  (73)  with  the  simulated  CADC  input  of  Mach  and  true 
airspeed.  Estimated  angle  of  attack  can  then  be  calculated. 

A  comparison  of  these  estimated  angles  of  attack  to  the 
documented  angles  of  attack  was  then  performed. 


\A 

* 


Table  III.  Static  Inflight  Estimator  Validation 

NT-33A  Aircraft 


DATA  POINT 

3 

4 

5 

6 

7 

8 

CW<DEG> 

.899 

-.899 

9.  402 

.802 

-.293 

2.498 

<WDEG) 

1.08 

-1 . 328 

9.198 

1.123 

-.712 

3.120 

Aerror 

(DEG) 

.  i76 

-.429 

-.204 

.321 

.414 

.622 

A* 


Table  III  depicts  the  err«,r  in  the  estimation  from  documented 
values.  With  Just  6  data  points,  the  average  error  over  the 
flight  range  is  .13  degrees  with  a  maximum  error  of  .622 
degrees  at  flight  condition  6.  However,  this  flight 
condition  is  the  only  data  point  at  40,000  feet,  resulting  in 
significant  least  square  error  in  the  regression  analysis, 
even  with  altitude  not  an  explicit  parameter  in  equation  (76). 
Throwing  out  that  data  point  results  in  an  average  angle  of 
attack  error  of  0.0337  degrees  over  the  remaining  range  of 
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1 


v-C}  flight  data  points.  This  regression  formula  is  vaiid  over 

■ «■» 

the  linear  portion  of  the  lift  curve.  It  can  be  extended  to 
the  stall  region  depending  on  the  extent  of  the  non-linear 
lift  region.  A  rapidly  stalling  wing-body  combination  could 
be  modeled  accurately  to  a* ,  while  a  flat  stall  region 
would  induce  angle  of  attack  errors  over  the  entire  range  of 
the  lift  curve. 


COMPUTER  TEST  PROGRAMS  Two  independent  methods  were  used  to 
verify  actual  angle  of  attack  and  sideslip  prior  to 
examination  of  the  Kalman  filter  operation  in  the  post-flight 
estimator.  These  programs,  called  LONG  and  LAT,  represent  a 
decoupled  system  for  modeling  the  longitudinal  and  lateral 
modes  of  response.  Input  into  LAT  and  LONG  consists  of 
flight  condition,  control  inputs,  initial  conditions  and 
stability  derivatives  of  the  body  axis  system.  Then  program 
then  forms  the  continuous  A  and  B  matrices  in  body  axis  and 
outputs  the  associated  eigenvalues  and  eigenvectors  for 
verification  with  actual  aircraft  data.  Forced  system 
response  to  step  inputs  are  calculated  using  a  Taylor  series 
expansion  to  second  order,  with  output  frequency  of  30 
samples  per  second  on  angle  of  attack  and  sideslip.  At  this 
point,  the  continuous  A  and  B  matrices  can  be  discretized  and 
input  into  the  discrete  Kalman  filter  programs,  DKF  and 
DKFLAT,  representing  the  longitudinal  and  lateral  estimators, 
respectively.  Using  the  same  step  control  inputs  as  in  LAT 
and  LONG,  the  a  and  (i  response  is  calculated  through  the 
transition  matrix  method.  These  responses  can  then  be 
compared,  and  if  formulated  corectly,  should  agree  to  an 
extremely  close  degree.  For  a  .3  degree  step  elevator  input, 
the  correlation  is  excellent,  as  shown  in  figure  9. 

Now,  DKF  and  DKFLAT  can  be  used  to  verify  the  ability  of 
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Figure  9.  Comparison  of  LONG  and  DKF  Derived  Angle  of  Attack 
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a  Kalman  estimator  to  recover  accurate  angle  of  attack  and 
sideslip  with  less  than  full  state,  noise-corrupted 
measurement.  Initially,  computer  generated  white  noise  was 
modified  to  a  intensity  equal  to  the  expected  error  in  the 
INS  measurement  signals.  This  Gaussian  noise  was  then  added 
to  the  exact  state  outputs,  and  these  states  became  noise 
corrupted  measurements  of  the  modeled  maneuver  for  which 
predicted  a  and  (i  traces  were  calculated.  The  Kalman 
filter  then  operates  on  the  available  state  measurements, 
estimating  the  remaining  states,  based  on  aircraft  model 
accuracy  and  measurement  uncertainty,  as  discussed  in  chapter 
IV.  Angle  of  attack  and  sideslip  were  then  calculated  from 
the  states  and  estimated  states.  Vhen  predicted  a  and  (i 
are  plotted  versus  Kalman  estimated  a  and  (i  ,  the  ability 
of  the  Kalman  filter  to  recover  these  parameters  becomes 
evident  and  quantifiable. 


POST-FLIGHT  VALIDATION  TEST  RUNS  The  post-flight  estimator 
uses  perturbations  from  trim  values  to  calculate  a  and  (i . 
As  such,  a  single  trim  flight  condition  was  used  as  a 
baseline  data  point.  This  condition  was  flight  condition  6 
of  Table  II,  using  the  NT-33A  aircraft  data. 


Table  IV.  Post-Flight  Validation  Test  Points 


DATA  POINT  1  2345 

INPUT  .  3°ELEV  2°ELEV  .3°RUD  2°RUD  2°AIL 


The  runs  produced  a  set  of  six  figures  (Figure  10 
through  Figure  13)  depicting  the  modeled  response  versus  the 


-  43  - 


Kalman  predicted  response. 


RESULTS  OF  VALIDATON  RUN  Initial  validations  runs  were  made 
with  the  INS  measurement  error  covariance  model  established 
in  Chapter  IV.  These  covariances  were  input  into  DICF  and 
DKFLAT  Q  and  R  matrices  to  determine  the  ability  of  the 
Kalivan  estimator  to  recover  angles  of  sideslip  and  attack 
with  noisy  INS  measurements.  For  these  initial  test  runs,  Q 
and  R  were  assumed  to  be  the  same,  reflecting  similar 
measurement  and  process  noises.  However,  the  accuracy  of  the 
INS,  as  Judged  by  the  covariances  that  were  calculated, 
provided  an  interesting  result.  The  Kalman  filter,  with  a  2 
degree  step  elevator  input,  calculated  the  exact  angle  of 
attack  as  the  model,  within  the  range  of  the  noise.  Figure 
10  depicts  this  result,  showing  the  estimated  angle  of  attack 
superimposed  on  the  modeled  angle  of  attack.  The  two  are 
identical  traces,  with  imperceptible  error.  This  is  hardly 
surprising,  considering  the  magnitude  of  the  measurement 
covariance  values  calculated  in  the  Post-Flight  Theory 
chapter.  The  noise  in  measurement  is  extremely  small.  These 
values  were  on  the  order  of  10” 8  radians  for  Euler  angles  and 
10 ~7  radxans/second  for  the  rates.  This  leads  to  the 

conclusion  that  indeed  angle  of  attack  and  sideslip  can  be 
accuately  estimated  with  less  than  full  state  measurement.  A 
properly  tuned  Kalman  estimator  can  make  up  for  less  than 
full-state,  noise-corrupted  measurement,  and  in  turn  allow 
for  accurate  estimates  of  a  and  /?,  with  an  equivalently 
accurate  INS. 

However,  a  worst  case  analysis  should  also  be  done  to 
allow  for  unknown  errors  and  noises  in  the  system  or  model. 
This  was  also  accomplished,  using  the  raw  accuracy  and  Jitter 
values  as  inputs  to  the  covariance  matrices  for  model  and 
measurement.  The  results  depict  an  accurate  system  as  well. 
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Figure  13.  Sideslip  Recovery  C.S  degree  Rudder) 
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though  these  covariances  are  3  orders  of  magnitude 
less  accurate  than  the  calculated  covariances.  Figures  11 
through  15  depict  these  worst  case  runs.  The  overall 
accuracy  of  the  estimations  is  at  the  worst,  . 5  degrees  in 
the  longitudinal  responses,  and  .4  degrees  in  the  lateral 
modes.  With  less  accurate  measurements,  as  this  scenario 
implies,  tuning  the  Kalman  filter  becomes  vastly  more 
important.  As  these  measurements  were  made  from  a  model, 
there  was  no  opportunity  to  actually  tune,  other  than  to  make 
a  good  initial  guess  of  model  accuracy.  This  best  guess  was 
the  input  of  measurement  covariances  into  the  model 
covariance  matrix.  While  this  is  a  good  starting  point  for 
covariance  analysis  of  a  physical  system  till,  a  real  world 
system  must  be  examined  to  tune  the  filter  accurately. 
However,  without  such  tuning,  both  modes  of  the  filter  are 
estimating  the  respective  angles  of  attack  and  sideslip  to 
about  1/2  degree.  A  best  guess  of  the  true  real  life 
performance  of  the  post-flight  estimator  probably  lies 
somewhere  between  the  caluclated  INS  covariance  analysis  and 
the  gross  input  of  INS  accuracy  values.  This  is  the  area 
that  must  be  examined  through  flight  test  with  actual 
ring  laser  gyro  equipment. 


VI.  FLIGHT  TEST 


TEST  PHILOSOPHY  The  purpose  of  the  flight  test  is  to  ensure 
proper  operation  of  the  inf] ight  angle  of  attack  routine 
using  actual  flight  data.  This  flight  data  collection  is 
critical  to  ensure  operation  with  noise-corrupted  input  data. 
The  noisy  nature  of  the  actual  INS  measurements,  coupled  with 
possible  unknown  or  assumed  effects,  will  test  the  ability  of 
the  algorithm  to  recover  the  necessary  data,  and  perform 
proper  computations  to  calculate  angle  of  attack.  Unlike 
simulated  data,  however,  no  true  angle  of  attack  is  known. 
Measured  angle  of  attack  from  calibrated  vanes,  normal 
aircraft  instrumentation,  and  computer  modeled  performance 
will  be  correlated  against  the  IMS  derived  a.  This 
V* »  information  will  provide  an  acceptable  measure  of  the 

accuracy  of  the  INS  derived  values  against  the  more 
conventional  approaches  in  obtaining  a. 

The  initial  riight  test  consisted  of  low  performance 
longitudinal  maneuvering  flight  data  tape  anlysis  only.  The 
reason  for  this  is  two-fold.  First,  it  will  demonstrate  the 
applicability  of  these  methods  for  a  determination  in  large, 
transport  type  aircraft.  As  these  aircraft  do  not  engage  in 
high-g  maneuvering  or  extreme  flight  attitudes,  the  basic 
assumptions  should  hold  throughout  the  nominal  flight  regime. 
The  applicability  of  these  determination  techniques  will 
be  demonstrated  for  large  aircraft  in  >oth  the  inflight  and 
flight  test  analysis  phases.  Secondly,  nominal  inflight 
accuracy  should  give  an  indication  of  the  proper  formulation 
of  the  estimator.  The  absence  of  high-g,  coupled  flight 
conditions  allows  a  straight  forward  evaluation  of  the 
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estimator.  The  robustness  of  the  estimator  in  maneuvering 
flight  will  be  discussed  in  Chapter  VII.  This  flight  test 
was  conducted  in  cooperation  with  NASA  and  in  conjuction  with 
a  NASA  propulsion  test  flight. 

TEST  AIRCRAFT  Initial  flight  test  was  accomplished  using  a 
NASA  F-15A  aircraft  manufactured  by  the  McDonnell  Douglas 
Corporation.  This  aircraft,  S/N  10281,  is  an  F-13A 
air-superiority  fighter  modified  for  digital  engine  and 
control  testing.  It  is  a  single  seat  aircraft  powered  by  two 
Pratt  and  Whitney  F-100  engines.  Flight  controls  consist 
of  twin  vertical  stabilizers  mounting  a  single  rudder  on 
each.  Lateral  control  is  effected  by  ailerons  on  the 
outboard  wing  surfaces,  aided  by  split  stabilators,  with 
pitch  controlled  by  symmetrical  stabilator  action.  The 
specific  aircraft  is  shown  inflight  in  Figure  16.  For 
further  information  on  the  basic  aircraft,  consult  T.O. 

IF— 15 A— 1  1171.  The  basic  aircraft  layout  and  critical 
dimensions  are  depicted  in  Appendix  E. 

INSTRUMENTATION  The  NASA  F-15A  Highly  Integrated 
Digital  Engine  Control  <HIDEC)  was  instrumented  with  a 
calibrated  yaw-angle  of  attack-pitot-static  CYAPS>  head. 

This  vane  provided  the  baseline  a  information  that  the 
inflight  estimator  was  evaluated  against.  In  addition,  the 
flight  test  boom  also  provided  pitot  pressure  for  calculation 
of  true  airspeed  within  the  central  air  data  computer  CCADC). 
This  aircraft  was  also  configured  with  production  angle  of 
attack  probes  which  were  used  as  a  secondary  comparison  of 
the  angle  of  attack  estimator.  The  external  instrumentation 
of  the  F-1SA  HIDEC  is  depicted  in  Figure  17. 

The  aircraft  was  fitted  with  a  USAF  standard  inertial 
navigation  system.  This  system,  a  Litton  Systems,  Inc. 
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ASN-109  inertial  navigation  set,  is  a  fully  self-contained 
dead  reckoning  navigation  system.  It  continuously  computes 
aircraft  position  by  double  integration  from  a  known  starting 
point.  Aircraft  ground  speed  and  attitude  are  interim 
computations  prior  to  position  computation.  This  INS 
consists  of  three  major  components.  The  first  is  the  actual 
inertial  measuring  unit  CIMJ)  which  houses  the  gyros  and 
accelerometers  to  sense  aircraft  motion.  The  second  is  the 
IMU  mount  which  provides  precision  mounting  and  alignment  of 
the  system  to  the  aircraft  body  axis.  The  third  component  is 
the  navigation  control  indicator  which  interfaces  the  INS  to 
the  central  computer  of  the  aircraft  and  also  allows  pilot 
control  of  the  functions  of  the  system.  The  INS  was  fully 
instrumented  and  a  list  of  signals  available  is  contained  in 
Appendix  F. 

Weight  was  available  through  production  fuel  sensors  on 
board  the  aircraft  which  measure  fuel  remaining  in  each  tank 
to  an  accuracy  of  200  lbs.  The  weight  of  the  aircraft  could 
then  be  easily  calculated,  knowing  basic  aircraft  weight, 
serviced  fluid  weight,  stores  weight,  and  the  changing  fuel 
weight. 

Mach  number  and  altitude  signals  were  obtained  from  the 
Sperry  AN/ASIC-6  air  data  computer,  along  with  true  airspeed. 
Air  density  was  calculated  through  the  standard  exponential 
atmosphere  equation  for  input  into  the  primary  estimator 
equation.  Overall,  no  signals  were  used  which  would  not  be 
obtainable  through  current,  INS-equipped  production  aircraft 
instrumentation. 

All  required  stability  derivatives  were  modeled  from 
flight  test  data  obtained  in  USAF  technical  reports.  C,  , 

w 

C  ,  center  of  gravity  motion  and  all  momenta  of  inertia  data 
mo 
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plots  are  included  in  Appendix  G.  In  addition,  all 
computer— derived  models  are  also  included  for  comparison 
purposes  in  the  same  appendix. 

DATA  REDUCTION  All  inflight  data  was  reduced  using  the  NASA 
ELEXI  computer  system.  Aircraft  data  telemetry  was  retrieved 
from  the  computer  for  the  specific  maneuvers  required.  This 
data  consisted  of  time  tagged  values  for  all  signals 
specified  in  the  estimator  flow  chart.  The  basic  estimator 
program,  as  implemented  in  FORTRAN  77,  was  altered  to  allow 
use  in  the  time  tagged,  sequential  mode  of  operation  on  the 
ELEXI  system.  This  consideration  was  important  in  that  the 
inflight  estimator  was  designed  for  real  time  operation,  and 
the  ELEXI  provided  that  capability  in  reducing  flight  data. 
Once  the  data  was  calculated  sequentially,  the  results  were 
plotted  and  compared  to  YAPS  boom  a  at  the  same  time  tag. 
FORTRAN  coding  of  the  inflight  a  estimator  for  use  on  the 
ELEXI  is  included  in  Appendix  H. 

TEST  METHODS  AND  CONDITIONS  The  optimum  flight  test 
technique  for  stable  longitudinal  flight  at  varying  angle  of 
attack  was  determined  to  be  the  level  acceleration.  In 
general,  the  aircraft  was  stabilized  on  conditions  in  the 
slow  speed  regime  with  engines  at  the  planned  military  or 
maximum  power  settings.  This  procedure  required  a  climbing 
ent^y  to  the  test  point.  The  aircraft  was  then  allowed  to 
accelerate  to  its  maximum  speed  while  maintaining  constant 
altitude  and  one-g  flight.  This  required  a  constant 
reduction  in  angle  of  attack  throughout  the  level 
acceleration  maneuver. 

Three  level  acceleration  test  points  were  planned  to 
evaluate  estimator  angle  of  attack.  Additionally,  Mach  and 


altitude  effects  were  investigated  using  military  and  maximum 
accelerations  at  three  different  altitudes.  Level 
acceleration  test  points  are  summarized  in  Table  V. 


TABLE  V. 

LEVEL  ACCELERATION  EVALUATION 


TEST  POINT 


ALTITUDE 

POWER 

20,000 

MIL 

10,000 

MIL 

40,000 

MAX 

MACH  BAND 


The  next  logical  step  in  the  buildup  process  to  evaluate 
the  inflight  a  estimator  is  to  introduce  abruptness  into  the 
estimation  process,  while  still  restricting  maneuvers  to  the 
longitudinal  modes  within  the  plane  of  gravity.  The 
wings— level,  constant  g  pitch-up  flight  test  technique  was 
considered  the  optimum  for  this  phase  of  the  flight  test. 

This  technique  required  the  aircraft  to  be  stabilized  at  a 
constant  aim  altitude  and  Mach  number.  The  aircraft  was 
abruptly  pitched  to  a  series  of  positive  and  negative 
constant  g  values,  much  like  a  roller  coaster.  This 
technique  was  accomplished  within  a  standard  2,000  foot  data 
band.  Three  test  points  were  evaluated  during  this  phase  and 
are  summarized  in  Table  VI. 

The  final  stage  of  the  quantitative  flight  test 
evaluation  of  the  angle  of  attack  estimator  is  examination  of 
out-of-plane  maneuvers.  The  purpose  of  this  phase  of  testing 
is  to  remove  the  gravity  vector  from  the  longitudinal  plane 
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TABLE  VI. 


ABRUPT  PITCH  EVALUATION 


TEST  POINT 

ALTITUDE 

MACH 

AIM  LOAD  FACTOR 

4 

40,000 

0.7 

•*•3.0 

5 

40,000 

0.7 

-0.3 

6 

40,000 

0.7 

-1.0 

of  the  aircraft  and  Judge  the  effect  of  load  factor  and 
banked  flight  of  the  estimator.  The  wind-up  turn  was  Judged 
to  be  the  best  flight  test  technique  for  this  flight  test 
goal.  In  this  maneuver,  the  aircraft  is  trimmed  at  a  given 
Mach  ,and  altitude.  The  aircraft  is  then  steadily  banked 
into  a  constant  Mach  turn  while  slowly  increasing  load  factor 
to  the  desired  end  point.  This  maneuver  takes  place  within  a 
2000  foot  data  band,  as  a  descent  is  required  to  maintain 
constant  Mach  at  a  trim  power  setting.  Two  test,  points  were 
identified  for  examination  of  the  effect  of  g  and  bank  angle 
on  the  estimator.  The  wind-up  turn  test  points  are  defined 
in  Table  VII. 

TABLE  VII. 

WIND-UP  TURN  EVALUATION 


TEST  POINT 

ALTITUDE 

MACH 

AIM  LOAD  FACTOR 

7 

20,000 

0.9 

+3.0 

8 

20 , 000 

0.  9 

+3.0 
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RESULTS  AND  ANALYSIS 


The  results  of  the  military  power  level  acceleration  in 
test  point  1  are  depicted  in  Figure  18.  This  acceleration 
was  from  0.5  to  0.9  Mach  at  20,000  FT.  The  data  is  presented 
as  two  angle  of  attack  traces.  The  first  is  AINF,  or  alpha 
infinity,  as  derived  fron  the  YAPS  boos.  The  second  is  AVB, 
or  alpha  wing-body,  representing  the  output  of  the  inflight  a 
estimator  based  on  equation  C17).  Also  accompanying  the  a 
traces  are  Mach  and  altitude  data  throughout  the  maneuver. 

This  data  will  accompany  all  level  acceleration  data  traces 
examined  in  this  section. 

There  are  two  areas  of  interest  in  Figure  18.  The  most 
obvious  is  the  initial  portions  of  the  a  data  traces,  where 
an  approximately  1.3  degree  noise  in  estimator  angle  of 
attack  is  apparent.  The  actual  level  acceleration  maneuver 
does  not  begin  until  17  seconds  into  the  data  trace.  This 
initial,  high  angle  of  attack  regime  is  the  climb  into  the 
maneuver  at  slow  flight.  This  slow  flight  is  characterized 
by  thrust  set  at  test  power,  in  this  case  military  setting. 
The  result  is  slow,  climbing  flight  in  moderate  buffet,  with 
some  internal  vibration  present.  The  difference  between  a 
sources  is  reasonably  constant  during  this  entry  into  the 
acceleration.  Initial  skepticism  of  the  estimator  would  give 
more  credence  to  the  YAPS  boom  a  than  the  estimator.  It 
should  be  realized  that  there  is  no  absolute  source  of  angle 
of  attack  in  this  test.  However,  the  a  estimator  would  be 
susceptible  to  airframe  buffet  and  vibration,  cluttering  the 
normal  load  factor  signal  at  these  low  speed  conditions.  It 
is  of  note,  though,  that  the  a  estimator  does  follow  the 
peaks  of  the  YAPS  boom  a  exactly,  remaining  within  1.3 
degrees  until  the  initiation  of  the  pushover  at  the  beginning 
of  the  level  acceleration,  occurring  at  13  seconds 
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into  the  trace.  The  YAPS  boon  appears  to  be  the  the  most 
accurate  source  during  this  slow  flight  phase. 

Once  the  level  acceleration  has  begun ,  the  two  angles  of 
attack  remain  within  0.5  degrees  of  each  other  throughout  the 
rest  of  the  maneuver.  For  the  first  half  of  the 
acceleration,  the  inflight  a  estimator  is  below  the  YAPS 
boom  a.  At  0.82  Mach,  the  traces  coincide,  with  estimated  a 
becominr  larger  than  boom  a  for  the  remainder  of  the  trace. 
They  do  stay  generally  within  0.5  degrees  during  this 
exchange.  Overall,  the  two  angle  of  attack  traces  coincide 
well,  with  the  exception  of  the  entry  into  the  maneuver, 
during  slow  flight  in  moderate  buffet. 

Figure  19  depicts  the  results  of  a  military  power  level 
acceleration  to  0.9  Mach  at  10,000  FT.  Spikes  in  this,  and 
subsequent  figures  indicate  data  dropout.  Again,  the  same 
two  phases  of  the  level  acceleration  are  notable  in  this 
figure.  During  the  slow  flight  entry  into  the  acceleration, 
angle  of  attack  traces  differ  by  approximately  2  degrees. 

The  inflight  a  estimator  follows  the  peaks  of  the  YAPS  boom  a 
exactly,  but  the  true  a  is  difficult  to  ascertain  for  this 
flight  regime.  However,  once  the  manuever  begins  at  20 
seconds  into  the  data  trace,  the  angles  of  attack  coincide 
well,  again  within  0.5  degrees.  This  is  the  attempted 
specification  to  which  the  estimator  was  designed  to  meet. 
Again,  the  crossover  of  a  traces  occurs  in  the  0.8  Mach 
regime,  as  was  noted  in  the  previous  level  acceleration. 

This  crossover  can  most  likely  be  attributed  to  the 
regression  of  the  lift  curve.  Recall  that  a  was  defined  as  a 
function  of  Mach  and  C,  .  At  0.82  Mach,  the  regression 
appears  almost  exact,  while  at  other  points,  there  is  some 
deviation  from  the  exact  lift  curve.  It  is  important  to  note 
at  this  point  also  that  no  divergence  occurs  as  the  higher 
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subsonic  Mach  numbers  are  approached,  indicating  a 
reasonably  good  fit  of  Mach  number  in  the  high  subsonic 
regime. 

The  maximum  power  level  acceleration  to  Mach  1.3  in 
test  point  3  also  showed  good  estitmator  correlation  to  YAPS 
boom  angle  of  attack  in  the  subsonic  regime.  Estimator 
tracking  can  be  seen  in  Figure  20,  including  Mach  and 
altitude  traces.  Pitch  maneuvers  tracked  well  through  almost 
17  degrees  boom  angle  of  attack.  Again,  the  aircraft  is  in 
slow  flight,  accompanied  by  moderate  buffet,  while  performing 
the  climbing  entry  to  the  level  acceleration.  This  is  the 
most  likely  cause  of  noisy  normal  accelerometer  output  at  the 
high  angles  of  attack.  Throughout  the  remainder  jf  the 
subsonic  portion  of  the  manuever,  both  a  traces  correspond 
nicely,  even  during  some  large  angle  of  attack  excursions.  A 
point  of  note  occurs  at  the  Jump  to  supersonic  flight  at  118 
seconds  into  the  trace.  At  this  point,  the  traces  begin  to 
diverge  at  a  rate  proportional  to  the  Mach  number.  The 
inflight  a  estimator  was  not  modeled  for  supersonic  flight 
and  this  could  be  the  simple  cause.  However,  the  traces  tend 
to  coincide  in  terms  of  deviations  from  a  steady  condition. 

In  other  words,  a  0.3  degree  jump  in  YAPS  boom  angle  of 
attack  is  matched  at  the  identical  time  segment  by  a  0.5 
degree  jump  in  estimator  angle  of  attack.  This  again 
indicates  proper  formulation  of  the  estimator,  but  any  YAPS 
boom  supersonic  errors  present  make  it  impossible  to  quantify 
estimator  errors  due  to  lack  of  supersonic  modeling.  In 
addition,  the  closeness  of  the  traces  indicate  at  this  point 
that  altitude  has  little  to  no  effect  on  the  model. 

Recall  that  the  approximation  was  made  to  eliminate  altitude 
from  the  model.  Initial  flight  test  shows  that  this  was  an 
acceptable  approximation  of  the  true  model. 
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Figure  20.  Maxiaua  Power  Level  Acceleration  Result* 
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Figure  21  depicts  the  results  of  the  abrupt  pitch 
evaluation.  The  upper  graph  shows  estimator/YAPS  angle  of 
attack  correlation,  while  the  lower  graph  shows  corresponding 
normal  load  factor.  This  figure  includes  all  abrupt  pitch 
test  points.  Of  immediate  note  is  the  closeness  with  which 
estimator  a  follows  boom  a  below  approximately  16  degrees 
angle  of  attack.  Negative  g  excursions  match  almost, 
identically.  The  lack  of  high  angle  of  attack  modeling  is 
the  cause  of  the  deviation  at  the  peak  of  the  high  g  points, 
as  was  noticed  in  the  initial  portions  of  the  level 
accelerations.  However,  the  close  correlation  of  the 
separate  angle  of  attack  sources  through  rapid  changes  in 
angle  of  attack  and  load  factor  do  support  the  basic  concept 
of  this  form  of  a  estimator.  It  can  indeed  accurately 
recover  angle  of  attack  with  at  least  0.3  degrees  of 
precision  in  upright,  purely  longitudinal  motion. 

The  wind  up  turn  evaluation  did  uncover  some  angle  of 
attack  deviations  in  the  estimator.  Figure  22  depicts  the  3 
g  wind  up  turn  results,  along  with  normal  load  factor 
achieved  in  the  maneuver.  In  this  figure,  there  is  some 
definite  deviation  during  the  sustained,  high  g  portion  of 
the  maneuver.  Although  the  traces  match  in  terms  of  peak 
locations,  they  differ  by  almost  a  degree  at  the  sustained  g 
point.  Although  the  YAPS  boom  a  does  not  provide  an 
absolute,  true  a,  it  should  be  the  weighted  preference. 
However,  the  boom  a  does  show  almost  0.5  degrees  worth  of 
noise  in  its  signal,  while  the  estimator  is  slighly  smoother. 
The  same  result  is  true  with  the  3  g  wind  un  turn  presented 
in  Figure  23.  YAPS  boom  angle  of  attack  is  consistently 
lower  than  estimator  a  at  the  higher  sustained  g  plateau, 
although  its  signal  is  much  cleaner  than  the  preceeding 
graph.  Again,  the  peaks  of  each  source  natch  well,  with  no 
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Figure  21.  Abrupt  Pitch  Maneuver  Results 
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lac  noticable  in  the  estimator.  The  lift  curve  used  as  the 
model  was  not  corrected  for  load  factor,  and  was  a  trimmed 
lift  curve.  This  could  account  for  the  deviations  at  the 
sustained  higher  load  factors,  and  indicate  a  requirement  for 
a  closer  wing-body  model  of  angle  of  attack.  Another 
possibility  is  error  in  out  of  plane  load  factor  calculation. 
This  could  be  the  result  of  actual  accelerometer  output  as 
opposed  to  theoretical  accelerations  about  the  center  of 
gravity.  The  result  is  an  alteration  in  equation  <283  to 
replace  the  cosOcos^  term  with  1.0  as  these  angles  are 
accounted  for  due  to  normal  accelerometer  bias  of  1  g.  This 
bias  is  included  in  all  normal  accelerometers  to  take  into 
account  the  gravitational  pull  of  the  earth.  In  straight  and 
level,  unaccelerated  flight,  the  normal  accelerometer  reads 
0.0  ft/sec2  acceleration  of  the  aircraft  center  of  gravity. 
However,  the  aircraft  is  indeed  under  32.2  ft/sec2  or  1  g 
acceleration  due  to  the  earth's  pull. 

Overall,  the  a  estimator  correlated  to  the  YAPS  boom  a 

well.  Under  most  conditions,  the  results  were  within  the 

specified  0.S  degree  deviation.  Where  the  deviations  were 

greater  than  that  value,  the  estimator  errors  were 

explainable  and  indicate  a  need  to  form  a  more  precise  C, 

Lw» 

model  than  a  linear  regression  on  two  independent  variables 
as  was  accomplished  for  this  research.  A  full,  maneuvering 
flight  demonstration  will  indicate  the  degree  to  which  the 
current  model  and  equations  are  adequate  for  high  g,  rolling 
flight  out  of  the  longitudinal  plane  of  motion. 
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VII.  MANEUVERING  FLIGHT  DEMONSTRATION 


PURPOSE  The  robustness  of  the  inflight  estimator  is 
evaluated  in  this  phase  of  the  flight  test  through  a  series 
of  highly  dynamic  maneuvers  in  varying  planes.  The  overall 
goal  is  to  highlight  weaknesses  in  the  inflight  estimator 
and  examine  regimes  of  flight  where  the  estimator 
assumptions,  as  currently  proposed,  break  down.  The  most 
likely  area  of  trouble  was  determined  to  be  out-of-plane,  or 
three-dimensional,  fighter  maneuvering  of  the  type  expected 
in  basic  air-to-air  or  air-to-ground  combat.  This  is 
therefore  the  emphasis  during  the  robustness  check  of  the 
inflight  a  estimator. 

SCOPE  The  purpose  of  this  portion  of  the  flight  test  was 
demonstration  only.  The  attempt  was  made  to  devise  a  single 
flight  test  technique  to  quickly  and  efficiently  demonstrate 
any  possible  area  of  weakness  in  the  angle  of  attack 
estimator.  In  other  words,  this  portion  of  the  investigation 
was  to  highlight  areas  to  troubleshoot  the  estimator 
algorithm  or  to  point  where  future  investigations  should  be 
directed. 

MANEUVER  The  robustness  check  was  only  accomplished  at  one 
flight  condition  due  to  constrained  flight  test  time.  The 
maneuver  that  was  developed  was  therefore  a  dynamic  one 
encompassing  all  expected  problem  areas  such  as  loaded 
rolls  and  longitudinal  pulls  out  of  the  local  horizontal 
plane.  The  modified  split-S  maneuver  was  performed  in 
conjunction  with  a  NASA  propulsion  test.  The  actual  NASA 
flight  test  card  is  included  in  Appendix  I. 

The  overall  robustness  maneuver  can  be  divided  into  four 


distinct  segments.  Initially,  the  aircraft  is  flown  in  a 
true  north  heading.  Once  established  on  conditions,  a  30 
degree  banked,  climbing  turn  at  2  g's  is  begun.  This  is 
indeed  a  climbing  turn,  as  a  level  30  degree  turn  requires 
only  1.2  g's.  Upon  stabilization  in  this  turn,  the  pilot 
then  rolls  inverted  in  the  same  direction  as  rolling  into  the 
2  g  turn  initially.  At  this  point,  the  pilot  then  begins  a 
sustained  4  g  pull  in  a  split-S  maneuver,  recovering  in  an 
upright,  wings  level  attitude.  The  robustness  maneuver  thus 
evaluates  a  climbing,  loaded  turn,  a  pure  roll  out  of  the 
local  horizontal  plane,  and  a  loaded  pull  with  the  gravity 
vector  constantly  moving  throughout  the  aircraft  axis  system. 
The  robustness  test  points  are  summarized  in  Table  VIII. 


TABLE  VIII. 
ROBUSTNESS  EVALUATION 


TEST  POINT 

ALTITUDE 

NACH 

ROLL  DIRECTION 

9 

20,000 

mwm 

RIGHT 

10 

20,000 

■ 

LEFT 

RESULTS  AND  ANALYSIS  The  results  of  the  first  robustness 
maneuver  are  depicted  in  Figure  24.  Pitch,  roll  and  yaw 
rates  are  presented  with  the  estimator  and  YAPS  boom  a 
traces.  The  initial  roll  into  the  maneuver  begins  at  2 
seconds  into  the  trace.  The  initial  difference  between  the 
higher  estimator  trace  and  the  YAPS  boom  trace  is 
approximately  0.7  degrees.  As  with  the  loaded  rolls 
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presented  in  the  Wind  Up  Turn  test,  the  difference  between 

traces  remains  reasonably  constant.  However,  as  the  roll  to 

inverted  attitude  begins  at  13  seconds  into  the  trace,  the 

curves  come  to  within  0.2  degrees.  The  4  g  split-S  shows 

very  close  correlation,  through  recovery  at  38  seconds.  Of 

note  during  this  phase  of  the  maneuver  is  the  noise  within 

the  a  estimator  signal.  While  no  more  than  approximately  0.8 

degrees,  it  disrupts  an  otherwise  close  match  under 

sustained  g,  inverted  flight.  The  second  robustness 

maneuver,  depicted  in  Figure  25,  shows  almost  the  exact  same 

results  for  the  opposite  direction  maneuver,  although  the 

loaded  turn  portion  is  more  distinct  in  this  plot.  Lack  of  a 

definite  bias  during  either  direction  of  the  maneuver 

indicates  that  sign  conventions  in  the  moment  correction 

equations  are  correct.  The  split-S  maneuver  in  the  second 

plot  also  shows  YAPS  boom  oscillations  of  up  to  1.5  degrees, 

while  the  estimator  is  smooth  in  relation.  Again,  the 

matching  is  excellent  during  the  recovery  phase  at  42 

seconds.  In  general,  the  only  deficient  area  of  the  inflight 

estimator  as  tested  is  the  bias  noticed  under  sustained  load 

factor.  Note  that  this  situation  did  not  occur  with  abrupt 

pitch  maneuvers.  This  deficiency,  on  the  order  of  0.3 

degrees  per  g  always  occurs  to  the  high  side.  Again,  a 

trimmed  lift  curve  at  1  g  was  used  as  the  model.  A  higher 

order  model  of  a  as  a  function  of  C  ,  altitude,  I!ach  and 

lvb 

load  factor  may  provide  the  key.  However,  correlation  during 
these  extreme  maneuvers  was  quite  acceptable,  considering  the 
multiple  changes  in  plane  and  velocity  vector  during  40 
seconds  of  robustness  evaluations.  The  maximum  difference 
was  2  degrees  as  observed  in  the  Wind  Up  Turn  test,  and  this 
occurred  under  approximately  5.5  g's.  In  addition,  it  is 
important  to  note  that  when  in  error  under  g,  the  a  estimator 
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was  always  higher  than  YAPS  boom  angles  of  attack. 


GENERAL  COMMENTS  The  ability  to  sit  in  the  NASA  control  room 
while  the  flight  was  in  progress  provided  a  unique 
opportunity  to  view  the  a  estimator  performance  during  all 
phases  of  flight.  While  it  is  difficult  to  quantify  all 
comments,  certain  qualitative  observations  can  be  made  from 
personal  engineering  notes  taken  during  the  flight. 

As  indicated  in  the  Flight  Test  chapter,  supersonic 
effects  were  evident.  In  general,  at  subsonic  speeds  during 
straight  and  level  flight,  the  u  estimator  was  well  within 
the  specified  0.5  degree  tolerance  as  compared  against  the 
YAPS  boom.  At  0.72  Mach,  the  estimator  was  within  0.3 
degrees  of  the  noseboom.  At  0.8  Mach,  the  estimator  tracked 
very  well  at  0.1  degrees  off.  And  x'inally  at  0.95  Mach,  the 
estimator  was  within  0.4  degrees  of  the  YAPS  boom.  At 
supersonic  speeds,  the  estimator  jumped  to  0.8  degrees  lower 
than  YAPS  values.  Again,  this  larger  difference  can  be 
attributed  to  the  purely  subsonic  modeling  of  the  lift  curve 
and  pitching  moment  of  the  F-13A.  However,  the  consistency 
of  the  low  estimator  values  indicates  that  the  inflight 
estimator  requires  some  "tuning”  to  better  approximate  the 
F-15A  aerodynamics. 

Two  secondary  aerodynamic  effects  were  observed  which 
deserve  note.  At  subsonic  speeds,  with  good  estimator 
correlation  to  YAPS  values  of  a  C within  .2  degrees)  , 
extension  of  the  massive  F-15A  speedbrake  caused  an  immediate 
Jump  to  0.3  degrees  difference  between  a  sources.  Two 
possible  explanations  are  readily  apparent.  First,  the  large 
speedbrake  alters  the  aerodynamic  characteristics  of  the 
wing,  invalidating  the  model  developed  for  clean 
configurations  only.  The  second  aerodynamic  effect  was  a  1.2 
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degree  estimator  to  YAPS  boom  difference  during  air  to  air 
refueling  at  320  knots  indicated  airspeed.  The  effect  here 
seems  to  be  caused  by  tanker  wake  effects  on  the  F-15A  local 
airflow,  impacting  both  the  YAPS  noseboom  and  the  local 
aerodynamics  as  modeled  by  the  a  estimator. 

Overall,  the  inflight  a  estimator  performance  was 
acceptable  as  qualitatively  evaluated  during  this  flight  and 
measured  against  YAPS  boom  values. 


VIII.  CONCLUSIONS  AND  RECOMMENDATIONS 


CENERAL  In  broad  terms,  the  objectives  of  this  thesis  were 
met  and  the  concept  of  angle  of  attack  and  sideslip 
estimators  using  standard  inertial  reference  platforms  is 
highly  feasible.  Of  the  two  estimator  variations,  the 
inflight  a  estimator  was  the  most  extensively  tested,  and  the 
most  widely  applicable.  It  served  to  demonstrate  that  the 
concept  of  an  angle  of  attack  estimator  which,  is  accurate  to 
0.5  degrees  is  not  only  possible,  but  available  for  real  time 
inflight  use,  with  current  generation  mechanical  inertial 
navigation  systems.  Specific  conclusions  and 
recommendations,  as  organized  by  thesis  objective,  follow: 

Objective  1.  The  linear  recursive  estimator  lends  itself 
well  to  use  as  an  angle  of  attack  and  sideslip  estimator. 

With  less  than  full  state  measurements,  the  estimators  can 
easily  determine  a  theoretical  value  for  that  missing  state, 
in  this  case  a  or  fi.  Developement  of  the  model  was  quite 
straightforward  and  required  no  extensive  mathematics  other 
than  formulation  and  discretization  of  the  A  and  B  matrices. 
Accuracy  achieved  by  this  system  is  difficult  to  measure. 
Using  the  noise  corrupted  computer  model,  a  0.25  degree 
accuracy  was  easy  to  achieve  with  ring  laser  gyro  accuracies 
and  variations  as  a  model.  However,  these  values  were  only  a 
first  guess  for  the  actual  measurement  covariance  matrix.  In 
addition,  model  uncertainties  were  likewise  only  first, 
guesses.  However,  the  system  can  be  tuned  with  a  computer  to 
provide  extraordinary  accuracies.  Currently,  the  system  as 
designed  is  accurate  for  determining  highly  accurate 
perturbed  angles  of  attack  and  sideslip  from  known  trimmed 
conditions.  The  system  should  be  expanded  to  incorporate 
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estimation  of  the  aerodynamic  angles  throughout  the  flight 
regime  and  independent  of  trimmed  conditions. 

R—l  DIRECT  FUTURE  STUDY  OF  a  AND  fi  LINEAR  RECURSIVE 
ESTIMATORS  TO  ELIMINATE  THE  REQUIREMENT  FOR  KNOWN  TRIMMED 
FLIGHT. 

Objective  2.  Both  angle  of  attack  and  sideslip  algorithms 
could  easily  be  validated  by  computer  simulation.  Indeed, 
the  computer  simulation  used  by  Gleason  and  modified  for  this 
effort  was  of  tremendous  value.  Given  discrete  A  and  B 
matrices,  a  full  state  simulation  could  be  run  simultaneously 
with  a  noise  corrupted  simulation  of  less  than  full  state. 

The  degree  of  noise  could  be  adjusted  as  required  to  closely 
model  the  actual  lateral  and  longitudinal  systems.  The 
ability  of  the  modeled  system  to  recover  a  and  f)  was  then 
successfully  demonstrated.  Again,  this  computer  simulation 
was  restricted  to  small  perturbations  about  a  trimmed 
condition.  In  addition,  the  solution  was  limited  to  only  one 
flight  condition  without  extensive  modeling  of  27  different 
parameters.  This  is  a  tremendous  limitation  to  its  use  in 
the  inflight  case  without  extensive  airborne  computational 
power. 

Objective  3.  An  inflight  angle  of  attack  estimator  was 
successfully  developed  for  use.  It  is  of  note  that  the  self 
imposed  requirement  for  minimum  calulations,  and  hence 
maximum  computational  speed,  did  not  restrict  the  accuracy  of 
the  estimator  as  tested.  Another  goal  under  this  objective 
was  also  reached.  All  signals  used  in  the  inflight  a 
estimator  were  from  standard  INS  or  onboard  data  sensors 
carried  by  almost  all  operations  military  aircraft.  In  other 
words,  apart  from  the  data  telemetry  systems,  no  special 
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flight  test  instrumentation  was  required  fcr  this  estimator. 

In  addition,  mechanical  INS  platforms  were  used,  allowing 

incorporation  of  this  estimator  in  current  generation 

aircraft  as  the  need  arises.  Accuracy  of  the  system  was 

highly  dependent  on  the  modeling  of  the  stability 

derivatives.  The  three  moments  of  inertia  proved  to  be 

secondary  effects,  not  requiring  extensive  mathematical 

models.  However,  C  was  a  critical  factor  and  was  used 

mo 

specifically  to  tune  the  system  to  YAPS  boon  angle  of  attack 
at  the  beginning  of  flight  test  data  evaluation.  This  step 
in  the  test  process  points  to  an  area  of  limitation.  An 
aircraft  still  must  undergo  some  flight  testing  with  a  YAPS 
boom  prior  to  using  an  a  estimator.  The  computer  program 
must  be  “calibrated"  to  the  aircraft,  at  least  in  the  initial 
flight  test  stages,  as  the  a  estimator  requires  historical 
data  to  model  the  lifting  system.  It  is  then  obvious  that 
the  most  critical  model  must  be  the  lift  curve,  with  a  as  a 
function  of  Mach  number  and  altitude.  The  most,  significant 
limitations  of  the  estimator  as  formulated  for  this  research 
was  the  lack  of  high  a  modeling  and  the  lack  of  a  supersonic 
capability.  This  limitation  was  based  solely  on  regression 
algorithms  available  to  formulate  this  estimator  version.  A 
more  powerful  regression  tool  would  allow  incorporation  of 
flight  regimes  that  were  not  modeled  by  the  current 
estimator. 

R-2  CONDUCT  SENSITIVITY  ANALYSIS  OF  C  ON  INFLIGHT  a 

mo 

ESTIMATOR  RESULTS. 

Objective  4.  The  inflight  a  estimator  was  evaluated  with 
flight  test  data,  and  robustness  examined  during  real-time 
flight  test.  Several  comments  can  be  made  as  a  result  of 
this  test.  The  primary  result  is  that  the  concept  is  indeed 


-  82  - 


feasable.  With  a  basic  multivariable,  linear  regression 
modeling  technique,  accurate  angle  or  attack  estimates  can  be 
made  in  all  attitudes  to  approximately  0.5  degrees,  as 
demonstrated  by  the  robustness  maneuvers.  Two  areas  of 
interest  need  to  be  highlighted.  Pitch  acceleration  was  not 
available  in  the  HIDEC  configuration.  If  it  were,  it  could 
be  expected  to  be  a  noisy  signal  due  to  the  algorithms  used 
in  its  calculation.  A  difference  of  pitch  rates  per  time 
interval  was  used  to  approximate  q.  This  proved  to  be  a 
satisfactory  approximation  of  the  term,  which  was  usually 
very  close  to  0.0.  The  overall  correction  ter*  in  the 
primary  estimator  equation  was  subsequently  small  when 
compared  to  the  other  terms,  and  could  actually  be  neglected 
with  only  limited  loss  of  accuracy.  Secondly,  the  rotation 
of  the  normal  acceleration  from  body  to  wind  axes 
demonstrated  the  difference  between  theoretical  equations  and 
reality.  The  rotation  equation  itself  was  based  of  c. g. 
accelerations  of  aQ„.  and  a„„.  In  straight  and  level 
unaccelerated  flight  in  the  wind  axis  system,  these  variables 
should  be  zero,  whether  the  aircraft  is  upright  or  inverted. 
However,  a  real  accelerometer  which  is  trimmed  to  read  0 
ft/sec2  aQ2  in  upright  unaccelerated  flight  will  read  64.4 
ft/sec2  in  inverted  unaccelerated  flight.  The  result  is  that 
the  correction  term  in  the  rotation  equation,  cosOcos^,  which 
represents  the  component  of  the  earth's  gravity  vector,  as  a 
correction  to  a  perfect  accelerometer,  is  taken  into  account 
by  the  real  accelerometer.  It  can  be  replaced  simply  by  the 
constant  one  g  acceleration  of  gravity,  as  the  rotational 
correction  is  automatically  applied  in  the  real  accelerometer 
readings. 

Objective  3.  The  robustness  maneuver  demonstrated  that  the 
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concept  of  a  simple,  efficient  angle  of  attack  estimator  was 
achievable.  The  estimator  was  accurate  to  within  the  0.5 
degree  desired  specification  with  two  exceptions.  First, 
under  sustained  higher  g  loadings,  the  estimator  accuracy  was 
degraded  proportional  to  the  loading.  This  indicates  that  a 
g  correction  term  needs  to  be  modeled  in  the  a  regression. 
Under  high  g,  the  estimator  was  always  high,  and  this  is  the 
more  favorable  of  the  possibilities.  Use  of  the  trim  C 

Si 

curves  could  be  the  cause  of  this,  and  simple  Modeling  under 

g  of  the  wing-body  C  should  suffice  to  correct  the 

La 

estimator  back  to  predicted  values. 

R-3  DIRECT  INFLIGHT  a  ESTIMATOR  EXPERIMENTATION  TO  LOADED 
FLIGHT  CONDTITION  RESEARCH. 

Finally,  in  the  calculation  of  a  known  singularity  in 

the  Euler  angle  rotations  was  reached  at  90  degrees  pitch 

angle.  This  situation  could  easily  be  rectified  by  revetting 

to  an  earlier  guess  of  a,  and  holding  that  guess  between  the 

80  to  90  degree  pitch  angle  phases  of  flight.  A  second 

solution  would  be  a  hold  register,  allowing  the  previous 

estimated  a  to  become  for  the  next  tine  segment.  This 

would  eliminate  the  need  for  Euler  angle  rotations  to  find 

ct  _ _ in  the  first  place,  and  seems  to  be  the  better 

ouess 

solution. 

In  conclusion,  the  area  of  angle  of  attack  and  sideslip 
estimation  is  an  exciting,  challenging  arena,  encompassing 
many  disciplines  of  Aeronautical  Engineering  and  Statistical 
Estimation.  Its  uses  are  bounded  only  by  inagination,  and 
its  possibilities  for  incorporation  into  current  aircraft  arc 
limited  only  by  desire. 


APPENDIX  A 


COMPUTER  PROGRAM  DKF 
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DKF 


The  program  DKF  uses  the  discretized  A  and  B  matrices 
of  the  aircraft  model  to  create  a  simulation  of  aircraft 
response  to  control  inputs.  This  response  is  ideal.  Noise 
inputs  for  the  model  and  the  measurement  are  added  and  a 
Kalman  estimator  then  utilizes  the  noise  corrupted  9  and  q 
values  as  simulated  measurements  to  study  the  ability  of  a 
linear  estimator  to  recover  the  remaining  two  states.  The 
subroutine  ALPHA  then  calculates  the  angle  of  attack  from 
the  estimated  values. 
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SSTOR  VIE:  2 

(]i|BMB<BWa»«aWiWBld|B«DKF .  FOR1*p*n*a*t,<p*0*n*p«p*n*1*a*n*n*0*11*1*0*^ 

C ****  DISCRETE  KALMAN  FILTER  PROGRAM 

c****  ANGLE  of  attack  estimator  version 
c****  BY  D. GLEASON  MODIFIED  BY  J.ZEIS 

6  AUG  86 


c*** 

X  -  A*XD  +  B*V  V 

N  NCO.QD 

c*« 

Z  -  H*X  -*■  V  V 

S  NCO , RD 

C*0*0*1 

XM  *  A*XP 

C*“*°« 

c«»« 

PM  »  A*PP*AT  B»Q*BT 

C*** 

XP  -  XM  +  KG*C  Z  -  H*XM  ! 

C*"** 

PP  »  C  I  -  KG*H  !  ’•PM 

C««* 

G  <nm«r»  mmnwinn* 

KG  «  FM*HT*t  H*PM*HT  +  R 

C*a*a«oMPMi«INPUT/OUTPUT  TAPES 

C  TAPE  1  =■  MATRIXC  INPUT  MATRICES  A,  B, H,Q,  R  AND  XPCOD,PPCODD 
C  TAPE  2  -  RANVC INPUT  PROCESS  NOISED 
C  TAPE  5  =«  RANVC  INPUT  MEASUREMENT  NOISED 
C  TAPE  7  «  MEASCOUTPUT  OF  KALMAN  FILTER  RESULTSD 
C  TAPE  8  -  SPC OUTPUT  FOR  PLOTTING  STATE  AND  STATE  EST. D 
C  TAPE  9  -  PPC OUTPUT  FOR  PLOTTING  COVARIANCES  AND  GAINSD 
C  TAPE  1 0* APC OUTPUT  FOR  ALPHA, INPUT  FOR  USPLOD 
C*a«omntaiB«OECLARATION  STATEMENTS 
C>*d**EXACT  SYSTEM 

REAL  XC4D , XDC4D 
REAL  ZC2D,UC2,600D 
REAL  AEC4 , 4D , BEC4 , 2D , HEC2 , 4D 
REAL  VEC2 , 600 D , SIGMAWC2D , QEC2 , 2D 
REAL  VEC2 , 600D , SIGMAVC2D , REC2 , 2D 
REAL  RANGEC4D 
O*«M0DEL  SYSTEM 

REAL  XMC4D , XPC4D 


REAL  AC  1,4) ,  BC  A . 2) ,  HC2 , 1) , QC2 , 2) , RC2, 2) 
REAL  ATCi.l) ,8TC2,4),HTC4,2) 

REAL  PMCi, l),PPC4,4),KGC4,2) 

REA^  VK1C4,4), VK2C4.1), VK3C4.4), VKACA.l) 
REAL  ALP AC 600 ,2),TIMC600) 

CHARACTER  TITLE*30 

INTEGER  IFLAGC1 ) , I , IY , N , M , INC, IOPT, IER 
C***Ko»**i*pgoGRAM  CONSTANTS 
IFLAGC 1 ) =60 
ITEST  "  IFLAGC 1)  +  1 
IPRINT  -  0 
DT=«1.  /30. 

ERRX=0 . 


ERRV=»0 . 


Cm*ww«INPMT  DISCRETE  ACNXN)  ,  BCNXM)  ,  HCLXN)  MATRICES 
CtfMamDfomMiNPUT  PROCESS  AND  MEASUREMENT  NOISE  COVARIANCE 

Ofrw>n».>».>.in»HATRICES  QCN,  M)  ,  RCL,  L) 

C-**INPUT  INITIAL  STATE  ESTIMATE.  XMCN)  ,  AND  ERROR  COV.PPCN.N). 
OPENC 1 ,  FILE3* ' MATRIX.  IN’) 

READC1,*)  N ,  M ,  L 
READC!,-*) 

READC 1 ,  -*)  CCAECI, J), J»1,N), I*1,N) 

RE ADC 1,*) 

READC  1 ,  -*)  CCACI, J), J«1,N), I-l ,N) 

READC1 , *> 

READC 1 , >*)  CCBECI, J) . J»1 , M) , 1-1 , N) 

READC  1 ,  >*) 

READC  1 ,  -*)  CCBCI.J),  J*1 , M) , 1-1 , N) 

RE ADC 1,*) 

READC 1 ,  -•)  CCHECI, J) , J-l , N) , 1-1 , L) 

RE ADC 1.*) 

READC 1,*)  CCHCI, J) . J»1 ,  N)  ,  1-1 ,  L) 

READC !,•) 

READC 1 ,  >*)  CCQECI, JD , J«1 ,  M)  ,  1-1 ,  M) 

RE ADC 1,*) 

READC  1 ,  -*)  CCQCI,  J)  ,  J-l , M) , 1-1 , M) 

RE ADC 1,#) 

READC1,-*)  CCRECI.J), J-l ,  L)  ,  1-1 ,  L) 

READC!,-*) 

READC1,-*)  CCRCI,  J)  ,  J-l ,  L)  ,  1-1 ,  L) 

READC!,*) 


v 

V 

V 

•J 

a 


* 

k' 

Yj 

$ 


a 


a 


a 


:4 
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READC  1  .  •*>  CXPCI)  ,  1*1  ,  N) 

READC 1 ,  *) 

READC i , *0  CCPPCI, J), J=1,N), 1=1, N) 


w*"*0*0"  INPUT  GAUSSI ANCO, 1)  PROCESS  NOISE 
0PENC2 . FILE= ' RANW. IN ' ) 

OPENC  7 , F ILE= ' MEAS . OUT ' ) 

OPENC  8 , FILE* ' SP . OUT ’ ) 

OPENC 9 , FILE* 'PP. OUT’) 

OPENCIO, FILE*’ AP. OUT’) 

READC 2  ,  "*)  NCYCLE 

VRITEC10,*)  NCYCLE 

VRITEC9,*)  NCYCLE 

VRITEC7 , *)  'NCYCLE  *’, NCYCLE 

VRITEC7 , *)  ’INPUT  PROCESS  NOISE’ 

READ  C2,*)  CCVECI.J),  1*1,  M),J*i,  NCYCLE) 

VRITEC7 , 105) C  C VEC I , J) , J*i , 10) , 1*1 , M) 

VRITEC7 ,  *) 

■"«"•»•  IN  PUT  GAUSSI ANCO ,  1 )  MEASUREMENT  NOISE 
WRITEC7 ,  *0  'INPUT  MEASUREMENT  NOISE' 

OPENC 5, FILE*' RANV.  IN') 

READC  5 ,  **) 

READ  C5,*)  C C VEC I, J),I»l,L),J*l, NCYCLE) 

VRITEC7, 103)  CCVECI, J) , J*1 , 10) , 1*1 , L)  VRITEC7 , ♦) 


■"•"•■•MODIFY  PROCESS  NOISE  TO  GAUSSIANCO ,  SIGMAV) 
■"•"•"•MODIFY  MEASUREMENT  NOISE  TO  GAUSSI  ANCO ,  SIGMAV) 
CALL  NOISECVE, VE, SIGMAW, SIGMAV, QE, RE, M, L, NCYCLE) 
VRITEC7,*)  ’MODIFIED  PROCESS  NOISE' 

WRITEC7 , 105)  CCVECI, J),J*1, 10), 1*1, M) 

VRITEC7 , <•) 

VRITEC7 , »)  'MODIFIED  MEASUREMENT  NOISE' 

VRITEC7, 105)  C C VEC I , J) , J»1 , 10) , 1*1 , L) 

WRITEC7,  >») 


INPUT  SEQUENCE  UCI,K) 


DO  20  K*l, NCYCLE 
DO  20  1*1, M 
UCI,K)-0. 

20  CONTINUE 
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DO  21  K*1.NCYCL£ 

UC1 .  JO  *0 . 00G7 
21  CONTINUE 

G>nln  irm  Iriln  |i  i<rHKn<n»>»i*n>n>n»Kn<p»Kn«n»iXB<nXBtn>^nXnxn>nxn<n»^n» 

C-»°«H*p«n«n>B)n>«^UTP»T  DISCRETE  A  B  AND  K  .MATRICES 

WRITEC7 ,  <0  \N  -  M  -  ' ,  M,  ’  L  -  ',L 

CALL  MPRINTC A£, N , N ,  'AE  MATRIX  ') 

CALL  ilPRINTC  A .  N ,  N ,  'A  MATRIX  '  > 

CALL  MPRINTC BE , N , M , 'BE  MATRIX  ’) 

CALL  MPRINTC B.N..M,  *B  MATRIX  ') 

CALL  MPRINTCHE, L, N, 'HE  MATRIX  ') 

CALL  MPRINTCH, u.N. 'H  MATRIX  ') 

C 

C>nnm»..«n»r'ir»ouTPUT  MEASUREMENT  AND  NOISE  COVARIANCE  MATRICES 
CALL  MPR1NTCQE,M,M, ’PROCESS  COV. -QE  ’> 

CALL  MPRINTCQ,M,M, 'PROCESS  COV. -Q  ’) 

CALL  MPRINTCRE,L,L, 'MEASUREMENT  COV. -RE  ') 

CALL  MPRINTCR.L.L, 'MEASUREMENT  COV. -R  ') 


C<«a«t*>**INITIALIZE  STATE  AND  ESTIMATE  VECTORS 
DO  50  1*1  ,N 
XC  I)*0. 

xdci)*o. 

50  CONTINUE 
DO  52  1*1, N 
XMC I  )*0 . 

52  CONTINUE 
DO  55  1*1 , N 
DO  55  J»1 , N 
PMC  I , J)»0 . 

55  CONTINUE 

CuninWFinmmw.mnoUTPUT  INITIAL  CONDITIONS  ON  X  &  PP 

VRITEC7,*)  ’mmmmmmmm *********  m mm  »»ww 

VRITEC7,*)  'CYCLE  -  0  TIME  -  0' 

CALL  VPRINTCXP.N, 'INITIAL  STATE  EST-XP') 
CALL  MPRINTCPP,N,N, 'INITIAL  CO VAR  EST-PP') 

C^bmbmdmbmniMAIN  LOOP  ON  K 
DO  1000  K*1 , NCYCLE 
TIME  -  !C*DT 

C hi iwwi >C ALCUL ATE  STATE  VECTOR 

CALL  STATECK,X,XD, VE,U, AE,BE,N,M,NCYCLE) 
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-CALCULATE  MEASUREMENT  VECTOR 

CALL  MEASURE  CK, X, Z, HE, VE, N, L, NCYCLE) 


C— CALCULATE  PREDICTED  ESTIMATE  &  PREDICTED  ERROR  COVARIANCE. 
CALL  XPREDCK, XM, XP, A, B, U, N, M, NCYCLE) 

CALL  COVARMC  A,  PP,  AT,  B,  Q,  BT,  N,  M,  PM,  VK1 ,  VK2  ,  VK3  ,  VK4) 

C 

C-CALCULATE  FILTER  GAIN, FILTER  ESTIMATE  &  FILTER  ERROR  COVARIANCE 
CALL  KGAINC PM , HT , H , R , N , L , KG , VK1 , VK2 , VK3 . VK4 ) 

CALL  XFILTERC XM , KG , Z , H , N , L , XP , VK1 , VK2 ,  VK3 ) 

CALL  CO VARPC  KG , H , PM , N , L , PP , VKi , VK2 , WK3 ) 

VRITEC9 , 105)  TIME, PPC1 , 1) , PPC1 , 2) , PPC2 , 2) 

VRITEC9, 105)  KG< 1,1), KGC  2 , 1 ) 


C— — —UPDATE  DELAY  VECTOR  XDCI) 

DO  120  I-1,N 
XDCI)«XCI) 

120  CONTINUE 

CALL  ALPH AC  X , XM , XP , N , K , ALFA , TIM , TIME) 

C— — PERFORM  ERROR  ANALYSIS 
C  CALCULATE  ERROR  INDICES 
EV»CXC2)-XPC2)) 

EX=*CXC1)-XP<D) 

ERRX-ERRX  +  ABSCXC1 )-XP< 1 ) ) 

ERRV-ERRV  ♦  ABSCXC2)-XPC2) ) 

C  WRITE  PLOT  VECTORS  TO  TAPE 

VRITEC8, 130)  TIME, XC2) , XPC2) , XC1) , XPC1 ) , EV, EX 
130  F0RMATC8F15. 5) 

C———— OUTPUT  FILTER  RESULTS 
IPRINT  *IPRINT  +  1 
IF  <. IPRINT  .EQ.  ITEST)  IPRINT  »  1 
IF  CIPRINT  .LT.  IFLAGC1))  GO  TO  1000 
VRITEC7 ,  — )  ' m  •  mw  . . . 

WRITEC7 ,  <•)  'CYCLE  «  ’,K,  '  TIME  -  ' ,  TIME 

CALL  VPRINTCX,N, 'STATE  VECTOR-X  ') 

CALL  VPRINTCZ,L, 'MEASUREMENT  VECTOR-Z') 

VRI TE  C  7 ,  -)  ' Mmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmmm. 

CX  CALL  VPRINTCXM.N, 'PREDICTED  EST. -XM  ') 

CX  CALL  MPRINTCPM,N,N, 'PREDICTED  COV. -PM  ') 
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CALL  VPRINTCXP.  N.  'FILTER  EST. -XP  ’> 

CALL  MPRINTCKG.N.L, 'KALMAN  GAIN  MATRIX  ’) 

CALL  MPRINTCPP, N. N. 'FILTER  COV. -PP  ’) 

WRITEC7  ,  »)  ’  wininirilninliirilnw  n» 

1GOO  CONTINUEC  wi.wm.il/vNGLE  OF  ATTACK  OUTPUTw"M»*»'»'»'»i»«.»ii"»»wnii 
VRITEC 1 0 ,  *)  ’TIME  EXACT  ALPHA 

CPREDICTED  ALPHA ' 

DO  551  K-t.NCYCLE 

VRITEC 10 ,549)TIMCK)  ,  ALPACK,  1  )  ,  ALP  AC  K .  2) 

549  F0RMATC2Y. F6 . 3, 16X, FI 2. 10.20X.F12. 10) 

551  CONTINUE 
C  OUTPUT  ERROR  INDICES 
ERRX=*ERRX/NCYCLE 
ERRV-ERRV/NCYCLE 

VRITEC  7 ,  *•)  'AVERAGE  POSITION  ERROR  -  ' ,  ERRX 
VRITEC7 , *)  'AVERAGE  VELOCITY  ERROR  -  ' , ERRV 
C<MMM°*M°«FORnAT  STATEMENTS 
100  F0RMATC2C  F15 . 5 . 2X) ) 

105  F0RMATC1"KF1  9  5.2X)) 

STOP 

END 

C 

CCCCCCCC1CCCCCCCC2CCCCCCCCC3CCCCCCCC4CCCCCCC5CCCCCCCC6CCCCCCCC7 
C  SUBROUTINE  VPRINT 

C  THIS  SUBROUTINE  PRINTS  OUT  A  VFCTOR  VITH  A  TITLE 
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 
SUBROUTINE  VPRINTC X , N , TITLE) 

REAL  XCN) 

CHARACTER  TITLEOO 
VRITEC 7 , 200) TITLE 
VRITEC 7,*)  X 
VRITEC7,*) 

200  F0RMATC25X , A20) 

RETURN 

END 

C 

CCCCCCCC1 CCCCCCCCCCCCCCC3CCCCCCCCC4CCCCCCCC5CCCCCCCC6CCCCCCCCC7 
C  SUBROUTINE  NPRINT 

C  THIS  SUBROUTINE  PRINTS  OUT  MATRIX  ACMXN)  VITH  TITLEC20  CHAR) 
CCCCCCCCCCCCCCCCCCCGCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 
SUBROUTINE  MPRINTCA,M,N, TITLE) 

REAL  ACM, N) 


"  n  " 


^"TC«ry_»  ».»«»■  W  V*  V  OS.”  W*-'  «-■  «-".  V"  f  MMC  V  V.  V  *.-  V'  V  ■<.'  «  "  * 
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CHARACTER  TITLE*30 
WRITEC7 , 2005TITLE 
DO  CO  1=1, M 

VRITEC7 , 300)  C AC  I , J) , J=1 , N> 

20  CONTINUE 
VRITEC7 ,  <0 
200  F0RMATC25X, A20) 

300  F0RMATC  8C2X , E9 . 3) ) 

RETURN 

END 

C 

CCCCCCCCC1CCCCCCCCCCCCC3CCCCCCCCC4CCCCCCCCC5CCCCCCCCC6CCCCCCCCC7 
C  SUBROUTINE  MADD 

C  THIS  SUBROUTINE  ADDS  TVO  MATRICES  CACMXN)  +  BCMXN)  »  CCMXN)) 
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 
SUBROUTINE  MADDCA, B, C, M, N)  REAL  ACM, N> , BCM, N> , CCM, N> 

DO  10  1=1, M 
DO  10  J=1,N 

CC I , J)=AC I , JD  +  BC I , J) 

10  CONTINUE 
RETURN 
END 


C 

C 

CCCCCCCCC1CCCCCCCCCCCCCC3CCCCCCCC4CCCCCCCC5CCCCCCCCC6CCCCCCCCC7 

C  SUBROUTINE  MTRANS 

C  THIS  PROGRAMS  TAKES  THE  TRANSPOSE  OF  MATRIX  ACMXN)  AND  RETURNS 

C  ATCNXM). 

CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 
SUBROUTINE  MTRANSC A, AT, M, N) 

REAL  ACM, N) , ATCN, M) 

DO  10  1=1, N 
DO  10  J=1,M 

ATCI, J)»ACJ,I) 

10  CONTINUE 
RETURN 
END 

C 

CCCCCCCCC1 CCCCCCCCCCCCCCC3CCCCCCCC4CCCCCCCC3CCCCCCCC6CCCCCCCC7 

C  SUBROUTINE  MMULT 

C  THIS  ROUTINE  MULTIPLIES  TVO  MATRICES. CACLXM)  X  BCMXN)  =  CCLXN)). 

C 

CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 
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SUBROUTINE  MMULTC A,  B ,  C,  L ,  M ,  ND 
REAL  ACL. MD , BCM, ND ,CCL, N> 

DO  10  1*1, L 
DO  10  J-l.N 
CCI, JD*0. 

DO  10  K-l.M 

CCI, JD-CCI, JD+  AC  I , KD-BCK, J) 

10  CONTINUE 
END 
C 

CCCCCCC1CCCCCCCCC2CCCCCCCCC3CCCCCCCC4CCCCCCCC5CCCCCCCC6CCCCCCCC7 
C  SUBROUTINE  NOISE 

CTHIS  SUBROUTINE  CHANGES  THE  GAUSSI ANCO , t D  PROCESS  AND  MEASUREMENT 
CNOISE  TO  GAUSSIANCO ,  SIGMAVD  ,  AND  GAUSSI  ANCO , SIGMAVD  RESPECTIVELY. 
CTHIS  CHANGE  IS  BASED  ON  THE  PROCESS  AND  MESUREMENT  NOISE 
CCO VARIANCE  MATRICES. 

CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 
SUBROUTINE  NOISEC  V, V, SIGMAV, SIGMAV,  Q ,  R,  M, L, NCYCLED 
REAL  VC M , NCYCLED , SIGMAVC MD ,QCM,MD 
REAL  VC L, NCYCLED, SIGMAVC LD,RCL,LD 
C— CALCULATE  NOISE  STANDARD  DEVIATION  FOR  PROCESS  NOISE 
DO  10  I«1,M 

SIGMAVC 1 D -SQRTC QC I,  I)) 

10  CONTINUE 

C— CALCULATE  NOISE  STANDARD  DEVIATION  FOR  MEASUREMENT  NOISE 
DO  20  I«1,L 

SIGMAVCID-SQRTCRCI, IDD  20  CONTINUE 
C— CREATE  GAUSSI ANCO , SIGMAVD  PROCESS  NOISE 
DO  30  1*1, M 
DC  30  J*1 , NCYCLE 

WCI, JD-VCI, JD-SIGMAVCID 
30  CONTINUE 

C— -CREATE  GAUSSIANCO , SIGMAVD  MEASUREMENT  NOISE 
DO  40  1*1, L 
DO  40  J-l, NCYCLE 

VC  I , JD -VC  I , JD -SIGMAVC I D 
40  CONTINUE 
RETURN 
END 

CCCCCCCCC1 CCCCCCCCC2CCCCCCCCC3CCCCCCCC4CCCCCCCC5CCCCCCCC6CCCCCCCC7 
C  SUBROUTINE  STATE 

C  THIS  SUBROUTINE  CALCULATES  THE  TRUE  STATE  VECTOR  AT  TIME  t 
C  XOCD-AE— XDCKD  +  BE-UCKD  ♦  BE-VCKD 
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CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 
SUBROUTINE  STATEC K, X , XD, V, U , AE , BE. N, M, NCYCLE) 

REAL  XCN) , XDCN) , WCM , NCYCLE) , UCM, NCYCLE) 

REAL  AECN , N) , BECN , M) 

DO  60  1=1 , N 
XCI)=0. 

DO  60  J=1 , N 

XCD-XCD+AECI,  J^XDCJ) 

60  CONTINUE 

DO  70  1=1, N 
DO  70  J=i,M 

XCD-XCD+BECI,  J)*VCJ,K:>-i-BECI,  J)*UCJ,ID 
70  CONTINUE 
RETURN 
END 

CCCCCCCCC1CCCCCCCCC2CCCCCCCCC3CCCCCCCC4CCCCCCCC5CCCCCCCC6CCCCCCCC7 
G  SUBROUTINE  MEASURE 

G  THIS  SUBROUTINE  CALCULATES  THE  TRUE  SYSTEM  MEASUREMENTS 
G  ZOO-HE^XCK)  +  VC  K) 

GCCCCCCGCCCCCCCGCCCCCCCCCCCGCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 
SUBROUTINE  MEASURECK, X, Z, HE, V, N, L, NCYCLE) 

REAL  XCN) , 2CL) , HECL, N) , VCL, NCYCLE) 

DO  80  1=1, L 
ZCI>=0. 

DO  80  J=1,N 
ZCD-ZCD+HECI,  J)*XCJ) 

80  CONTINUE 
DO  90  1=1, L 
ZCI)=ZCI)+VCI,K) 

90  CONTINUE 
RETURN 
END 

CCCCCCCCC1CCCCCCCCC2CCCCCCCCC3CCCCCCCC4CCCCCCCC5CCCCCCCC6CCCCCCC7 
G  SUBROUTINE  XPREDICT 

GTHIS  SUBROUTINE  PREDICTS  THE  STATE  VECTOR  PRIOR  TO  THE  MEASUREMENT 
C  XMC  K) = A*XP<  K~1 )  +  B*U 

CCCCGCCCCCCGCCCCCCCCCCCGCCCCCGCGCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 
SUBROUTINE  XPREDOC , XM , XP , A, B , U , N , M , NCYCLE) 

REAL  XMCN) , XPCN) , ACN, N) 

DO  10  1=1, N 
XMCI5-0. 

DO  10  J=1 , N 

XMCD-XMCI)  +  ACI,  JD«XPCJ) 


10  CONTINUE 
DO  20  1*1 ,  N 
DO  20  J-l ,M 

XMCD-XMCD+BCI ,  J5*UC  J,  K) 

20  CONTINUE 
RETURN 
END 
C 
C 

CCCCCCCCC1CCCCCCCCC2CCCCCCCCCCCCCCCC4CCCCCCCC5CCCCCCCC6CCCCCCCC? 
C  SUBROUTINE  COVARM 

CTHIS  SUBROUTINE  CALCULATES  THE  ERROR  COVARIANCE  MATRIX  PRIOR  TO 
CTHE  MEASUREMENT 

C  PMCK)*A*PPCK-1  )*AT  +  B*Q*BT 

CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 
SUBROUT I NE  CO VARM C  A , PP , AT , B , Q , BT , N , M , PM ,  PP1 ,  PP2 ,  PP3 , PP4 ) 
REAL  ACN , ND , PPCN , N),PMCN,N),BCN,M),QCM,M> 

REAL  ATC  N ,  N)  ,  BTC  M ,  N ) 

REAL  PP1 CN , N) , PP2CN , ND , PP3CN , M>, PP4CN, N) 

C  PPl-A^PPCK-1) 

C  PP2*A**PPCK-1)>«AT-PP1*AT 

C  PP3*B>K? 

C  PP4«B*Q*BT*PP3*BT 

CALL  MTRANSC A, AT, N, N5 
CALL  MMULTC A, PP, PP1 ,  N, N, N) 

CALL  MMULT C  PP1 , AT , PP2 ,  N ,  N ,  N  > 

CALL  MTRANSC B, BT, N, M5 
CALL  MMULTCB, Q, PP3 , N , M, M) 

CALL  MMULTC  PP3 , BT , PP4 ,  N ,  M .  N> 

CALL  MADDC PP2 , PP4 , PM ,  N ,  ND 
RETURN 
END 
C 

CCCCCC1CCCCCCCCC2CCCCCCCCC3CCCCCCCC4CCCCCCCC5CCCCCCCC6CCCCCCCC7 
C  SUBROUTINE  KGAlN 

C  THIS  SUBROUTINE  CALCULATES  THE  KALMAN  GAIN. 

C  KGCKD«PM*HT'*CH*PM'«HT  -  R!-l 

CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 
SUBROUTINE  KGAINC  PM , HT , H ,  R ,  N ,  L ,  KG ,  K1 , K2 , K3 , K4) 

REAL  PMCN, N) , HCL, N> , HTCN, L) , KGCN, L) , RCL, L) 

REAL  K1CN, L) , K2CL, L) , K3CL, L) , K4CL, L) 

REAL  VKC1305 
C  Ki*PM"*HT 
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G  K2=H*PM*HT*H*fCl 

G  JC3=H*PM*HT+R*IC2+R 

C  KA=CH*PM*HT+R!  INVERSE=K3  INVERSE 

C  KG*P*HT*CH*PM*HT+R! INVERSE*  K1*K1 

GALL  MTRANSCH, HT,  L, N) 

CALL  MMULTC  PM ,  HT ,  K1 ,  N ,  N ,  L) 

GALL  MMULTC  H ,  IC1 ,  IC2  ,  L ,  N ,  L) 

GALL  MADDOC2 , R,  K3 ,  L ,  L3 

CALL  GMINVCL , L ,  K3  ,  K4 , 0 , 0  ,  L) 

CALL  MMULTC K1 ,  K4  ,  KG,  N ,  L,  L} 

RETURN 

END 

G 

CCCCCCC1CCCCCCCCC2CCCCCCCCC3CCCCCCCC4CCCCCCCC5CCCCCCCC6CCCCCCCC7 
G  SUBROUTINE  GOVARP 

C  THIS  SUBROUTINE  CALCULATES  THE  ERROR  COVARIANCE  AFTER  THE 
C  MEASUREMENT  HAS  BEEN  MADE. 

C  PPCK)»[  II— KGCK)>*H!  ’•PM 

CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCGCCCCCCCCCCGCCCCCCCCCCCCCC 
SUBROUTINE  COVARPC  ICG ,  H ,  PM ,  N ,  L ,  PP ,  II ,  PP1 ,  PP2) 

REAL  KGCN, L3 , HCL, N) , PMCN, N) , PPCN, N) 

REAL  IICN, N) , PP1CN, N5 , PP2CN, N) 

G  PP1*  -KG*H 

C  PP2*IT— KG—H-II+PPl 

C  CREATE  IDENTITY  MATRIX  IICNXN3 
DO  10  1*1 , N 
DO  10  J*1,N 
IICI, J)*0. 

IFCI.EQ.JD  IICI,J3-1.0 
10  CONTINUE 

CALL  MMULTCKG, H, PP1 ,  N,  L, N) 

C  NEGATE  PP1  MATRIX 
DO  20  1*1, N 
DO  20  J-1,N 

PP1 C I , J}*  -PP1CI,J) 

20  CONTINUE 

CALL  MADDC II, PP1 , PP2 ,  N,  N) 

CALL  MMULTC PP2 , PM , PP, N,  N, N) 

RETURN 

END 

C 

CCCCCCC1CCCCCCCCC2CCCCCCCCG3CCCCCCCC4CCCCCCCCSCCCCCCCC6CCCCCCCC7 
C  SUBROUTINE  XFILTER 
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C  THIS  SUBROUTINE  CALCULATES  THE  STATE  VECTOR  ESTIMATE  AFTER 
C  THE  MEASUREMENT  OCCURS. 

C  XPCK)=XMCK)  +  KGCZ-H*XMCK)  ! 

CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 
SUBROUTINE  XFILTERC XM , KG , Z , H , N , L , XP , HXM , RZ , KGR) 

REAL  XMCN)  , KGCN, L) , ZCL) , HCL, N) , XPCN) 

REAL  HXMC  L, 1 ) , RZCL, 1 ) , KGRCN , 1 ) 

C  HXM-H*XM 

RZ=Z-H«XM*Z-HXM*RESIDUALS 
KGR*KG*C  Z-H*XM !  *KG*R 
CALL  MMULTCH, XM, HXM, L, N, 1 ) 

C  NEGATE  HXM 

DO  10  1*1,  L 

HXMC  I,  1)  —  HXMC  1,1) 

10  CONTINUE 

CALL  M ADDC  Z , HXM , RZ , L , 1 ) 

CALL  MMULTCKG, RZ, KGR, N,  L, L) 

CALL  MADOCXM ,  KGR,  XP,  N ,  1 ) 

RETURN  END 


* 

*  SUBROUTINE  ALPHA 

»  THIS  ROUTINE  CALCULATES  THE  EXACT, MEASURED  AND 
■*  PREDICTED  ANGLE  OF  ATTACK,  AND  STORES  THE  VALUES 

*  IN  ARRAYS.  VTRIM  IS  THE  TRIM  VELOCITY,  AND  ALFT 

*  IS  THE  TRIM  ANGLE  OF  ATTACK. 


SUBROUTINE  ALPHAC X , XM , XP , N , K, ALPA, TIM , TIME) 
REAL  XCN),XMCN), XPCN) , TIMC600) 

REAL  ALPAC600 , 2) 

VTRIM-570.0 
ALFT*. 014 

UO»VTRIM*COSCALFT> 

VO* VTRIM  OSINCALFT) 

ARG1  *C  VO+XC  2  )  )  XC  UO+XC  1 )  ) 

ARG3-C  VO+XPC  2) ) /C  UO+XPC 1 ) ) 

C  ALPHA  CALCULATION 

ALPACK, 1 )*ATANC ARG1 ) 

C  ALPHA  PREDICTED  CALCULATION 
ALP  AC  K ,  2  )  -  AT  AN  C  ARG3  ) 

TIMCK)-TIME 

RETURN 


END 

SINCLUDE:  'GMINV.  FOR' 
^INCLUDE:  ’ VADD. FOR' 
^INCLUDE:  'DOT. FOR' 
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SUBROUTINE  BETA 


The  program  DKFLAT  is  identical  to  DtCF  with  the 
exception  of  the  subroutine  BETA  instead  of  ALPHA  as  in  DJCF. 
Input  and  output  is  change  to  accomodate  the  different 
variables,  but  all  logic  is  the  same.  BETA  will  calculate 
the  modeled  and  predicted  angle  of  sideslip  from  the  lateral 
equations  of  motion  as  described  in  the  program  and  text. 


SUBROUTINE  BET 

THIS  ROUTINE  FINDS  ACTUAL  AND 

PREDICTED  ANGLE  OF  SIDESLIP  AND  STORES  THE  VALUES 
IN  ARRAYS.  VTRIM  IS  THE  TRIM  VELOCITY,  AND  ALFT 
IS  THE  TRIM  ANGLE  OF  ATTACK.  * 


SUBROUTINE  BET<X,  XM, XP, N, K, BETA, TIM, TIME) 

REAL  XCN) , XMCN) , XP<N) , TIM<600) 

REAL  BETAC600 , 2) 

VTRIM-570.0 

ALFT=.014  J 

UO-VTRIM>*COSCALFT) 

VO-VTRIM  1*SIN<  ALFT) 

ARGl-XCl) /VTRIM 
ARG2-XP<1)  /VTRIM 
C  BETA  CALCULATION 

BETAC K, 1 ) »/SINC ARG1 ) 

C  BETA  PREDICTED  CALCULATION 
EETACK, 2)»ASIN(ARG2) 

TIMCK5-TIME 
RETURN 
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APPENDIX  C 


COMPUTER  PROGRAM  LONG 


LONG 


The  program  LONG  provides  an  independent  source  to 
compute  modeled  angle  of  attack  to  verify  DKF  derived  angle 
of  attack.  This  program  takes  the  continuous  stability  and 
control  matrices  and  calculates  through  a  Taylor  series 
expansion  the  response  of  the  system  to  step  control  inputs. 


••LONG.  FORCIOOO  DATA  POINT  VERSION)' 
wiwngiaciOO  SECONDS  MAX)  . . . »"** 


•^AERONAUTICAL  ENGINEERING  PROGRAMS 
•■BY:  J.E.  ZEIS  VERSION  5C9  APRIL  86)« 


«w^»«TIME  RESPONSE  OF  AIRCRAFT  < FORCED  AND  UNFORCED)’ 
IN  THE  LONGITUDINAL  MODES 

1.  30DY  AXIS  SYSTEM 
2-  FLIGHT  PATH  ANGLE  IS  SMALL 

3 .  XVDOT=0 . 0 

4 .  XQ=0 . 0 


mu ii  n  ir.mwmm.DATA  STORAGE  INSTRUCTIONS . . . ■■ .  ■  * 

1.  DATA  FOR  DIFFERENT  AIRCRAFT  IS  STORED  IN  SEPERATE  FILES 

i.e.  F-4C.DAT  OR  A-4D.DAT 

2.  THIS  PROGRAM  USES  PLANE. DAT  AS  ITS  INPUT  FILE.  COPY  AN 
AIRCRAFT  FILE  AS  "PLANE. DAT".  THEIR  FORMATS  ARE  EXACTLY 
THE  SAME. 

3.  THE  OUTPUT  FILE  IS  "LONG. PRN"  TO  SEE  THE  OUTPUT  ON 
YOUR  TERMINAL,  "TYPE  LONG. PRN". 

4.  TO  CHANGE  THE  TIME  OVER  WHICH  RESPONSE  IS  COMPUTED,  JUST 
CHANGE  TIME  IN  FILE  "PLANE. DAT". 


REAL  MU, MVDOT, MW, MQ, ME, MT 

DIMENSION  AC4 , 4) , BC4 , 2) , EAT<4 ,4), EC 4, 4), GCA, 4), AINVC4 , 4) 
DIMENSION  P<4 , 4) , PPC4 , 4) , ASUBC4 , 4) 

DIMENSION  U<2005) , VC2005) , THETC2005) . QC2005) , CONC4)  , 

C ALP AC 2005) 

DIMENSION  EATIC4 , 4) 

DIMENSION  EVRC4) , EVIC4) , VECRC4 , 4) , VECIC4 , 4) 

DIMENSION  INDICC4) 

DATA  TAU/O. 0333333/ 
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300  OPENC 2, FILE- 'PLANE. DAT') 

310  RE ADC  2 . 320) V, ALF, T J, XU , XV, XE , XT 

320  F0RMATC/,/,/,9X,F6.  1 , 28X,  F7 . 3  ,  /,  5X ,  F9. 6  f  /,  5X ,  F9.  6  ,  ✓  ,  5X , 

CF7 . 4 , 10X , F7 . 1 , /, 10X ,  F9 . 6 , /) 

330  READC2 , 340)ZQ, ZU, ZVDOT, ZV, ZE,  ZT 

310  F0RMATC5X , F7 . 3 , /, 5X , F7 . 4 , ✓, 9X , F9 . 6 ,  /,  5X , F7 . 1 , X, 10X , F7 . 2 , ✓ 

C10Xf F7. 2./) 

350  RE ADC 2 , 360 ) MU , M VDOT , M  V , MQ , ME , MT 

360  F0RnATC5X,F9.6,/,9X,F9.6,/'\5X,F9.6,/'  5X,  F5. 2,  /,  10X.  F7.  2  .  ✓ 

C10X,F9. 6,/) 

370  RE  ADC  2 , 380  )  UC 1 )  ,  VC  1 )  ,  THETC 1 )  ,  QC  1 )  ,  ELEV,  THROT  „  TI 

380  FORMATC/, /, 3X , F7 . 2 , ✓  , 5X , F7 . 2 ,  ✓,  9X , F7 . 3 , /, 5X, F7. 3 , /, /, /, 6X , 

CF7 . 4 , /, 6X , F7 . 1 , /, ✓, 33X , F5 . 1 ) 

VELOCITY  AND  TIME  *w°t°M^oMn<pwn»<wp» 

400  VO»V*SINCALF) 

410  uo-v*coscalf) 

420  IT-TI*30+1 

500  AC1 , 1 )— XU 

510  AC1 , 2)— XV 

520  AC1 , 3)*— 32 .194 

530  AC  1,4)— VO 

540  AC2 , 1 )— ZU/C1 . 0— ZVDOT) 

550  AC 2, 2)-ZV/Cl . 0— ZVDOT) 

560  AC2 , 3)-0 . 0 

570  AC2,4)«CV+ZQj/C1.  O-ZVDOT) 

580  AC3 , 1 )— 0 . 0 

590  AC  3 , 2 ) — 0 . 0 

600  AC 3, 3)— 0. 0 

610  AC3,4)-1.0 

620  AC 4, 1)»ZU*MVD0T/<1.  O-ZVDOT)  ♦  MU 

630  AC4,2)«ZV*MVD0T/C1.  O-ZVDOT)  +  MV 

640  AC 4 , 3)— 0 . 0 

650  AC4 , 4)»C V+ZQ)*MVD0T/C1 .  O-ZVDOT)  MQ 

Cnnwtur.inmiii.i.wiBwinwiB-MATRIX  CALCULATIONS1**1*1  *"«r  *’•"»  mmv  i . mm  . . . 

660  BC1,1)-XE 

670  BC1,2)-XT 

680  BC 2,  D-ZE/Cl.  O-ZVDOT) 

690  BC2,2)-ZT/C1. O-ZVDOT) 

700  BC3,1)»0.0 

710  BC3,2)-0.0 

720  3C 4, 1)-ZE*MVD0T/C1. O-ZVDOT)  ♦ME 

730  3C  4 , 2 ) -ZT*M VDOT/C 1 . O-ZVDOT)  +MT 


rf*  .■Vw'TIL’V  "s.  «  V.*  msi.*r 


.■>. 

W 


# 


C1*1*n*n*n»TntP<ptP<o4o<n<i  IDENTITY  MATRIX  INPUT*"*******  r  »  l>iln|r|nln|n>i|,nn|n».«n)n«ntn<n«n«. 

450  DATA  E/1.,.0,.0,.0,.0,1.,.0,.0,.0,.0,1.,.0,.0,.0,.0,1./ 

C)MUwa°w°w°w°«OPEN  OUTPUT  FILE  LONG. 

100  OPENC 4  ,  FILE=  ' LONG.  PRN  ’ ) 

C+tmxcmxrnrnnumnnnicuttSERIES  EXPANSION  FOR  UNFORCED  RESPONSE*********** 
C**»n«iilri«n«n|n«n<n»n|o«a«EXPC  AT)  =  I+AT+1/2CASQ])CTSQ)  +  .  .  .  ******************* 
1000  DO  1040  1=1,4 

1010  DO  1030  J=1 , 4 

1020  GCI,J)=ACI, 1 5*AC1 , J5+AC I , 2)*AC2 ,  JD+AC 1 , 3}*AC3  ,  J)+AC 1,4)* 
CAC  4 , J  3 

1030  CONTINUE 

1031  J=i 

1040  CONTINUE 

1050  DO  1110  1=1,4 

1060  DO  1100  J=1 , 4 

1070  EATC I ,  J  )=EC  I ,  J  )+TAU*AC  I ,  J)  +  .  5*TAU*TAU*GC  I ,  J) 

1100  CONTINUE 

1101  J=1 

1110  CONTINUE 

C****************CALCULATION  OF  “A  INVERSE“^>">"W>n*’«"«n»"«»«"»''»"W'»>'i'  >»*■— 

3000  NR=4 

3010  NC=4 

3020  MT=0 

3Q30  MR=0 

3040  NC0L=43045  DO  3049  1=1,4 

3046  DO  3049  J=1 , 4 

3047  aSUBCI, J)=ACI, J) 

3049  CONTINUE 

3050  CALL  GMINVCNR, NC, ASUD, AINV, MR, MT,  NCOL) 

C*****************CALCULATION  OF  "EAT  —  I ’* *********************** 

3200  E ATI Cl , 1)=EATC1 ,  ID— 1 . 0 

3210  E ATI Cl , 2)=EATC1 , 23 

3220  EATIC1 , 3)=EATC1 ,3) 

3230  EATIC1 , 4)=EATC1 , 4) 

3240  EATIC2, 1)=EATC2, 1) 

3250  EATIC2 , 2)=EATC2 , 2>-l . 0 

3260  EATI C  2,3) =EATC  2 , 3 ) 

3270  EATIC2 , 4)=EAT<2 , 4) 

3280  EATIC3 , 1 )=EATC3 , 1 ) 

3290  EATIC3 , 2)=EATC3 , 2) 

3300  EATIC3 , 3)»EATC3 , 3)— 1 . 0 

3310  EATIC3 , 4)=EATC3 ,4) 

3320  EATIC4 , 1)=EATC4 , 1 ) 
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3330  EATIC4,23-EATC4,23 

yj  3340  EATIC4 , 33*»EATC4, 3) 

3350  EATI C  4 , 4  3  =“EAT C  4 , 4  3  — 1 . 0 

OmwiwmwiwwnmninirniHULTIPLY  AINV  TIMES  EATI***  »■»■»'►  ^t^,,*************** 

3400  DO  3450  1-1 , 4 

3410  DO  3440  J-l , 4 

3420  PC  I ,  J3*AINVC  1 , 1 3 —EATI Cl ,  J3+AINVC  1 , 23—EATIC2  ,  J3+AINVC  1,3) 

C— EATIC3 ,  J3+AINVC 1 , 43—EATIC4 ,  J3 

3440  CONTINUE 

3441  J-l 

3450  CONTINUE 

Omin>n>r.»rl«rl>nln»Hn|.)n|..|n>iMULTIPLY  P  TIMES  B  MATRI X1>l«"«"»r  »">■>"»" *********  ****** ** 

3510  DO  3600  I=*l ,  4 

3520  DO  3590  J»1 , 2 

3530  PPCI,  J3=PCI,13*EC1,  J3+PCI.23-BC2,  J3+PC I , 33-BC3 ,  J3+PCI,43* 

CBC4 , J3 

3590  CONTINUE 

3591  J-l 

3600  CONTINUE 

Oih»"hi.»iiih»i.  --------FINAL  MULTIPLICATION  FOR  AINV-CEAT  -  I3-B-UCT: 

4000  C0NC13-PPC1,13-ELEV+PPC1,23-THR0T 

4010  CONC23»PPC2,13-ELEV+PPC2,23-THROT 

4020  C0NC33»PPC3,13*ELEV+PPC3,23»THR0T 

4030  CONC43-PPC4,13-ELEV+PPC4,23-THROT 

C  lull  I  ilnln—li——— —frill— STATIC  OUTPUT'!'  WM  -H—— Mi'—li  Wmmi.»,|.in»> 

8000  VRITEC4, 80103 

8010  F0RMATC5X , ’FORCED  AIRCRAFT  RESPONSE  TO  ELEVATOR^THROTTLE 

C  INPUT ’3 

8020  VRITEC4, 80303 

8030  F0RMATCX,8X, 'INITIAL  CONDITIONS  ARE  AS  FOLLOWS: '3 

8040  VR1TEC4 , 8050  3  UC1 3 , VC 1 3 , THETC 1 3 , QC 1 3 

8050  F0RMATC4X, ’DELTA  U- ' , F8 . 3 , 4X, 'DELTA  V-' , F8. 3 , 4X, 'THETA- ' , 
CF8.  4  ,  / ,  4X,  '  Q*»  ’ ,  F8 . 43 

8051  VRITEC4 , 80523ELEV, THROT 

8052  FORMATC/, 'ELEVATOR  INPUT  IS: ',F6.3,4X, 'THROTTLE  CHANGE: ', 

C  F6.13 

8060  VRITEC4 , 80703 

8070  FORMATC/, 6 X, 'STATE-SPACE  SYSTEM  OF  THE  FORM  XCD0T3-AX  ♦  BU'3 

8080  VRITEC4, 809038090  FORMATCX, 23X, ' "A”  MATRIX ' 3 

8100  DO  8130  1-1,4 

8110  VRITEC  4, 81203  AC  1, 13,  AC  I, 23,  AC  I, 33,  AC  1, 43 
a}.  8120  FCRMATC3X,  F10. 6 , 4X, F10 . 6 , 4X, F10. 6 , 4X, F12. 63 

^  8130  CONTINUE 
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8140  VRITEC4 ,3150) 

8150  FORMATC/, 23X,  ’  “A  INVERSE”  MATRIX’) 

8160  DO  3190  1=1,4 

8170  VRITEC  4 , 81 80) AINVC 1,1), AINVC 1,2), AINVC 1,3), AINVC 1 , 4) 

81 80  FORMATC  3X , FI 0 . 6 , 4X , F10 . 6 , 4X , FI 5 . 3 , 4X , FiO . 6) 

8190  CONTINUE 

8200  VRITEC4, 8210) 

8210  FORMATC/, 16X,  ”’B"  MATRIX’) 

8220  DO  8250  1=1,4 

8230  VRITEC 4 ,8240)B(I,1) ,BCI,2) 

8240  FORMATC 3X, FIO. 6,4X,F10. 6) 

8250  CONTINUE 

Cl*|,*l+l<0*‘*1*rl*n*rl*n*n*n*n*~l*rl*rltrE  I  GEN  V  ALU  E  AND  EIGENVECTOR  CALL***»»i>Mn<nininMi 
8500  N=4 

8510  NM=4 

8520  CALL  EIGENCN, NM, A, EVR, EVI , VECR, VECI , INDIC) 

8530  WRITEC4 , 8540) 

8540  FORMATC/, /,8X, ’EIGEN  ANALYSIS’) 

8550  VRITEC  4,8560) 

8560  FORMATC/, 1 OX, ’EIGENVALUES’) 

8570  DO  860C  1=1,4 

8580  VRITEC 4  3590)EVRC I) , EVIC I) 

3590  FORMATC 5X, E13. 6 , 2X, ’+’,E13.6, ’i’) 

8600  CONTINUE 

8610  VRITEC 4 , 8620) 

8620  FORMATC/, 1 OX, ’EIGENVECTORS’) 

8630  DO  8700  1=1,4 

8640  VRITEC 4, 8641) I 

8641  FORMATC/, 5X, ’VECTOR  *’ , 13) 

8642  DO  8699  J-l , 4 

8643  VRITEC 4 , 8644) VECRC J, I) , VECIC J, I) 

8644  FORMATC 2X, Ell .  4 , IX, ’  +  ’ ,  Ell .  4, ’ i ’ ) 

8699  CONTINUE 

8700  CONTINUE 

C'*»»‘*l<i|<n>n|n>nB.>nln>nlnlrnn»Hn<n|ipREQENCY  AND  DAMPING  RATIO»*°«°*°*^«*^^^ 

8710  0MEG1=SQRTCEVRC1)**2  ♦  EVIC1)*«2) 

8720  0MEG2=SQRTCEVRC3)**2  +  EVIC3)**2) 

8730  IF  C0MEG1  . LT. 0MEG2)  GOTO  8800 

8740  0MEGP=0MEG2 

8750  D AMPP= ABSC  EVRC  3 ) /OMEGP) 

8760  0MEGS=0MEG1 

8770  DAMPS-ABSCEVRCD/OMEGS) 

8780  GOTO  8850 


8800  0MEGP-0MEG1 

8810  DAMPP-ABSC£VRC1)/0MEGP) 

8820  0MEGS=0MEG2 

8830  DAMPS-ABSC  EVR(3)/0MEGS) 

8850  VRITEC 1 , 8860)0MEGS, DAMPS 

8860  F0RMATC/.5X, 'OMEGAC SHORT  PERIOD)- ' , El 1 . 4, IX, 

G 'SHORT  PERIOD  DAMPING  » ' ,  El 1 .  A  ) 

8870  VRITECA , 8880)OMEGP, DAMPP 

8880  F0RMATC5X, 'OMEGACPHUGOID)- ' , Eli . A, AX, 'PHUGOID  DAMPING 

CE11. A,/) 

OUTPUT  HEADER***** *** 

9000  VRITECA , 9010) 

9010  F0RMATC/.8X, 'AIRCRAFT  RESPONSE  (STABILITY  AXIS  SYSTEM)') 

9020  VRITEC A, 9030) 

9030  F0RMATC5X,  'T’,7X,  'U',11X,  'V',8X,  'THETA', 8X,  'Q',8X,  'ALPHA') 

C******************^^^  TIME  RESPONSE********"*!  * »"»**■ «"*"»****** 

1500  N *1 

1510  T-N*TAU 

1520  M-N+l 

C*i"»ii.i.ii.«ni,mii>.iiiiXCNfl)-EXPCAT)*XCN)  immmtmmmmmmmmmmmmmmm  . . . 

2000  UCM)*EATC1 , 1  )*U(N)+EAT(1 , 2)*VCN)*EATC1 , 3)*THETCN)+ 

C  EATC1 ,  A)*K|CN) 

2010  VCM)»EATC2,1)*UCN)+EATC2,2)*VCN)+EATC2,3)*THETCN)  + 

G  EATC2 , A)*QCN) 

2020  Q(M)-EATCA,l)*U(N)+EATCA,2)*VCN)-*-EATCA,3)*THETCN)* 

C  EATCA , A)*0CN) 

2030  THETCM)  =*EATC  3 , 1 )*U(N)+EAT(  3 , 2)*VCN)+EATC3 , 3)*THETCN)  + 

C  EATC3 , A)*QCN) 

2100  CONTINUE 

C+m*n*m*m*»o***mFO RGED  TIME  RESPONSE*********  *********  . . 

C**"*I****P**I**"*B*|XC  N+l )—XCN+l )  >  CONTROL  VALUE** »*wi***mwi^ww».».ii.»<inii 

2200  UCM)-UCM)+C0NC1) 

2210  VCM)-VCM)+C0NC2) 

2220  THETCM)-THET(M)+C0NC3) 

2230  QCM)=*QCM)+CONCA) 

22A0  ARG»CV0*VCM))XCU0+UCM)) 

2250  ALPA(M)-ATANCARG) 

C****  ****** . "INCREMENTAL  TIME  RESPONSE  OUTPUT  ****** . 

9030  VRITEC A, 9060)  T, UCM) , VC M) , THETCM) , QCM) , ALPACM) 

9060  F0RMATC2X,  F5 . 2 ,  AX,  F8 .  3 , 3X,  F7 . 3 , 3X,  F7.  A ,  AX,  F6.  A ,  AX,  F7.  5) 

9070  N-N+l 

9080  IF  CN.LT.IT)  GOTO  1510 

9090  CONTINUE 
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9100  END 
SINCLUDE:  'EIGEN. FOR' 
SINCLUDE:  ' GMINV. FOR’ 
SINCLUDE:  ' VADD. FOR’ 
SINCLUDE:  ’DOT. FOR' 
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"jL 


.  * . 


APPENDIX  D 

COMPUTER  PROGRAM  LaT 
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LAT  is  the  lateral  response  equivalent  of  the  program 
LONG.  Again,  the  continuous  stability  and  control  matrices 
are  input,  along  with  the  prscribed  step  control  input  to 
find  the  modeled  response  of  the  aircraft  in  the  lateral 
modes.  This  is  used  to  verify  the  operation  of  the  lateral 
estimator  program,  DKFLAT. 


»LAT.  F0RC1000  DATA  POINT  VERSION)' 


^"•AERONAUTICAL  ENGINEERING  PROGR  - 
►BY:  J.E.  ZEIS  VERSION  1C8  APRIL  86)> 


RESPONSE  OF  AIRCRAFT  C FORCED  AND  UNFORCED)' 
IN  THE  LATERAL  MODES 


ii.mmwr.ini.il  ASSUMPTIONS  *mtm 

1.  BODY  AXIS  SYSTEM 

2.  ASSUME  SMALL  FLIGHT  PATH  ANGLE. 

3.  ASSUME  YAV  RATE  IS  SMALL. 

4.  YCV  DOT)«L'<V  DOT)«N’CV  D0T)-0 


*"*T»n«DATA  STORAGE  INSTRUCTIONS!  .  i  im  >i  i  . . mmim 

1.  DATA  FOR  DIFFERENT  AIRCRAFT  IS  STORED  IN  SEPERATE  FILES 

i.e.  F-4C.LAT  OR  A-4D. LAT 

2.  THIS  PROGRAM  USES  PLANE. DAT  AS  ITS  INPUT  FILE.  COPY  AN 
AIRCRAFT  FILE  AS  "PLANE. DAT" .  THEIR  FORMATS  ARE  EXACTLY 
THE  SAME. 

3.  THE  OUTPUT  FILE  IS  "LATERAL. PRN"  TO  SEE  THE  OUTPUT  ON 
YOUR  TERMINAL,  "TYPE  LATERAL. PRN". 

4.  TO  CHANGE  THE  TIME  OVER  WHICH  RESPONSE  IS  COMPUTED,  JUST 
CHANGE  THE  TIME  IN  THE  LAST  LINE  OF  THE  INPUT  PROGRAM 
PLANE. DAT. 


REAL  LB , LR , LP , LVDOT , LA , LRD , NB , N VDOT , NP , NR,  N A ,  NRD 
DIMENSION  AC 4 , 4 ) , BC4 , 2) , EATC4 ,4) , EC4 , 4) , GC4, 4) , AINVC4, 4) 
DIMENSION  PPPC4 , 4) , PPC4 , 4) 

DIMENSION  VC2005) , PC2005) , PHIC2005) ,  RC2005) , C0NC4) , 
CBETAC2005) 

DIMENSION  EATIC4 , 4) , ASUBC4 , 4) 

DIMENSION  EVRC4) , EVIC4) , VECRC 4,4) , VECIC4,4) 

DIMENSION  INDICC4) 

DATA  TAU/O. 033333/ 

ii.iinmininmmmmirn  read  INPUT  LOGIC  . mt  i  wwmw  mi  i  mm 
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•V .« 


* 


u 

L* 


300  0PENC2 . FILE* ' PLANE.  DAT ’ D 

310  RE ADC 2 , 320) VEL, ALF, YP, YB, YR, YA, YRD 

320  FORMATC/, /,/,9X,F6.i,28X,F7.3,/,/,40X,F9.6,/, 40X, F9 . 4 , /, 4uX, 

CF7 .  4 ,  / , 15X,F7.4,/,45X,F9.6,/) 

330  READC2 , 340 )LR, LP, LVDOT, LB,  LA, LRD 

340  FORMATC 41 X ,  F7. 3 , /, 41X, F7 . 4 , /, 45X, F9 . 6 ,  /, 41X, F7.  2,/,  46X, 

CF7. 2 , /,  46X,  F7.  2 ,  /) 

350  READC2 , 360 )NB, NVDOT, NP, NR,  NA, NRD 

360  FORMATC 41 X , F9 . 6 , /, 45X, F9 . 6 , /, 41X, F9 . 6, /, 41X, F5.  2,  /, 46X, 

CF7. 2,/,46X,F9. 6,/) 

370  READC2 , 380) VC1 ) , PCI ) , PHIC1 ) , RC1 ) , AIL, RUD, TI 

380  FORMATC/, /, 40X, F7 . 2 , /, 40X, F7. 2 , /, 42X , F7. 3, /, 40X,  F7. 3 ,  /,  /, 

C/, 41 X, F7. 4 , /, 41X, F7 . 4 , /, /,  33X, F5 . 1 ) 
CiMr<nMiMr.|r^»<nW^If)ITIAL  VELOCITY  AND  TIME*"**"**^****^^ 

400  VO=»VEL*SINCALF) 

410  UO=VEL*COSCALF) 

420  IT=»TI*30+1 

C%i«nm»i»nn>Mir'«ii»»«i»»«mi»mi«ii«A--MATRIX  CALCULATIONS’>l,«"»»*l>»,»— * wx^twi  mi  mi  imin» 

500  AC1,1)-YB/U0 

510  AC1 , 2)=*W0 

520  AC  1 , 3)=*—32.  194 

530  AC1,4)»-U0 

540  AC  2 , 1 )  =*LB/UO 

550  AC2,2)*LP 

560  AC2,3)»0.0 

570  AC2,4)~LR 

580  AC3,1)«0.0 

590  AC3 , 2)=*1 . 0 

600  AC3 , 3)=0 . 0 

610  AC3,4)«0.0 

620  AC4,1)*NB/U0 

630  AC4,2)*NP 

640  AC4,3)*0.0 

650  AC4,4)*NR 

C»nn»Hn»iiDiaiuin>»w»»mnn»B-MATRIX  CALCULATIONS  H'Mimwinni«*»tB«a«atBir mmmmmmtatatafm 


660 

BC1,1)»YA 

670 

BC1,2)«YRD 

680 

BC2,1)»LA 

690 

BC2,2)*LRD 

700 

BC3,1)»0.0 

710 

BC3 , 2)=0 . 0 

720 

BC4,1)-NA 
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730 


450 
C* 

100 

C«**« 

c* 

1000 
1010 
1020 

1030 

1031 
1040 
1050 
1060 
1070 
1100 
1101 
1110 
C4BMM 

3000 
3010 
3020 
3030 
3040 

3046 

3047 

3049 

3050 
C* 

3200 

3210 

3220 

3230 

3240 

3250 

3260 

3270 

3280 

3290 

3300 

3310 


3C4,2)*NRD 

wp^o^ideNTITY  MATRIX  INPUT* 

DATA  E/1 .  ,.0,.0,.0,.0,1.,.0,.0,.0,.0,1.,.0,.0,.0,.0,1./ 
)PEN  OUTPUT  FILE  RESPOND. OUT* 

0PENC4,FILE-’LATERAL. PRN') 


“•"•’•“•‘•“•“•SERIES  EXPANSION  FOR  UNFORCED  RESPONSE* 


*EXPCAT)-I+AT+1/2<AS0)CTSQ)  +  .  .  .  *********** 

DO  1040  1-1,4 
DO  1030  J-1,4 

GCI, J>«A<I.1)*AC1 .  J)+ACI,2:>*AC2,  J>+AC 1 , 3)*AC3 ,  J)+ACI,4)* 
CAC4, J) 

CONTINUE 

J-l 

CONTINUE 
DO  1110  1-1,4 
DO  1100  J-1,4 

EATCI ,  J)-EC  I ,  J}-*-TAU*AC  I ,  J)  +  .  5*TAU*TAU*GCI,  J) 

CONTINUE 

J-l 

CONTINUE 

■  ■‘■■'CALCULATION  OF  "A  I N VERSE "  . mm i—  wwi 


NRV-4 

NC—4 

flR-0 

MT-0 

NC0L-43045  DO  3049  1-1 , 4 
DO  3049  J-l , 4 
ASUBC I , J)*AC I , J3 
CONTINUE 

CALL  GMINVCNRW, NC, ASU3, AlNV, MR, MT, NC0L) 
:ALCULATI0N  OF  "EAT  -  Vmmmum 
EATIC1 , D-EATCl , 1)-1 . 0 
E ATI Cl , 2>— EATCI , 25 
EATIC1 , 33-EATC1 ,3) 

EATIC1 , 4D-EATC1 , 4) 

EATIC2 , 1 }»EATC2 , 1 ) 

EATIC2 , 25-EATC2 , 25-1.0 
EATIC2 , 35-EATC2 , 3) 

EATIC2, 43-EATC2, 4) 

EATIC3 ,  D-EATC3 , 1) 

EATIC3 , 23— EATC3 , 25 
EATIC3 , 3D-EATC3 , 35-1.0 
EATIC3 , 4> -EATC3 , 4) 


\ 
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3320  EATIC4. 1)=£AT04,1) 

3330  EATIC4 , 2) *EAT<4 , 2) 

3340  EATIC4 , 3)=EAT<  4,3) 

33S0  EATI<4 , 4)»EAT(4 , 4)— 1 . 0 

AINV  TIMES  EAT  I  **m*m*^m*'n*m*nmm*m*'^^ 

3400  DO  3450  1*1,4 

3410  DO  3440  J*1 , 4 

3420  PPPCI , J)*AINVCI , 1 )*EATI<1 , J)+AINVC I , 2)*EATI<2 , J)+AINV< I , 3) 

C—E/TIC3 , J)+AINV< I , 4)*EATIC4 , J) 

3440  CONTINUE 

3441  J*1 

3450  CONTINUE 

Cxo«u«mn»j»*««n«»i«n«m*mflULTIPLY  ppp  TIMES  B  M ATR I X itolomatammtatMamatBlatMotalBtam 

3510  DO  3600  1*1,4 

3520  DO  3590  J*l,2 

3530  PPCI,  J)=PPP<I,l)*Bd,  J)+PPP<I,2)*B<:2,  J)+PPPC I , 3)*B<3 ,  J) 

C+PPPCI,4)*BC4, J) 

3590  CONTINUE 

3591  J»1 

3600  CONTINUE 

MULTIPLICATION  FOR  AINV*CEAT  -  I)*B*UCT)  — 
4000  CONd  )»PPd ,  l)*AIL+PPd  ,  2)*RUD 

4010  C0N<2)=PPC2,1)-AIL+PPC2,2)*RUD 

4020  C0N<3)»PPC3,1)-AIL+PP<3,2)>*RUD 

4030  C0N<4)*PPC4,1)*AIL+PP<:4,2)*RUD 

———————STATIC  OUTPUT— 

8000  VRITEC4 , 8010) 

8010  F0RMATC5X, 'FORCED  AIRCRAFT  RESPONSE  TO  AILERON/RUDDER 

C  INPUT’) 

8020  VRITEC4 , 8030) 

8030  F0RMAT</,8X, 'INITIAL  CONDITIONS  ARE  AS  FOLLOWS:') 

8040  VRITEC4, 3050) VC1), PCI ),PHId),Rd) 

8050  F0RMATC4X, 'DELTA  V* ' , F8 . 3 , ✓  , 4X, ' DELTA  P- ' , F8.  3 , 4X, 'PHI*' , 
CF8. 4 ,  /,  4X,  '  R* ' ,  FO .  4 ) 

8051  VRITEC4, 8052) AIL, RUD 

8052  FORMATC/, 'AILERON  INPUT  IS: ' ,  F6 . 3 , 4X, ' RUDDER  INPUT  IS:’, 
CF6. 3) 

8060  VRITE<4 , 8070) 

8070  F0RMAT</,6X, 'STATE-SPACE  SYSTEM  OF  THE  FORM  XCDOT)»AX  +  BU') 
8080  VRITEC4 , 8090) 

8090  FORMATC/, 23X,  ’ ”A"  MATRIX’) 

8100  DO  8130  1-1,4  \ 

3110  VRITEC4 , 8120) ACI , 1) , AC I , 2) , AC I , 3) , A< I , 4) 
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8120  FORMATC3X,  F10 .  6  ,  IX ,  F10 . 6 , IX ,  F10 . 6  ,  IX ,  F12.  6) 

8130  CONTINUE 

8110  WRITEC1 .8150) 

8150  FCRMAT</,23X, '"A  INVERSE"  MATRIX ' ) 

8160  DO  8190  1-1,1 

8170  VRITEC1 , 8180) AINVC 1,1), AINVC 1,2), AINVC 1,3),  AINVC  1,4) 

8180  F0RMATC3X ,F14. 6 , IX, FI 1. 6 , IX, FI 5. 3, IX, FI 1.6) 

8190  CONTINUE 

8200  VRITEC1 , 8210) 

8210  F0RMATC/.16X, '"B"  MATRIX') 

8220  DO  8250  1-1,1 

8230  WRITEC1 , 8210)BC I , 1 ) , B< I , 2) 

8240  F0RMATC3X, F10.6,1X,F10.6) 

8250  CONTINUE 

C*<r-»«n«nin>l>nan|iHn»— — —m—QUTPUT  HEADER**0**0* -*11*" — Ki-nmr- K—nniwtmiH  m>m>wn<nw- 

9000  VRITEC1 , 9010) 

9010  F0RMATC/,8X, 'AIRCRAFT  RESPONSE  (BODY  AXIS  SYSTEM)') 

9020  VRITEC 1,9030) 

9030  FORMATC5X,  'T',7X,  'V'.llX, ’P',8X,  'PHI  ’ , 8X, 'R' , 8X, 'BETA' ) 

Giminriinnir.ini  ■mm,  i  immum  ni.UN FORCED  TIME  RESPONSE  ■  . .  —mu—  in m —m— 

1500  N— 1 

I  1510  T-N-TAU 

1520  M-N+l 

2000  VCM)«EATC1 ,  1  )>*VCN)*EATC1 ,  2)*PCN)-*-EATCl ,  3)>*PHICN)+EATC1 ,1) 

C*RvN) 

2010  PCM)=EATC2,^)‘*VCN)+EAT<2,2)*PCN)+EATC2,3)*‘PHICN)+EATC2,1) 

C-RCN) 

2020  RCM)— EATC1 , 1  )’*VCN)-»,EATC1 , 2)’*PCN)+EATC1 , 3)>*‘PHICN)-*,EATC1 , 1) 

C*-RCN) 

2030  PHICM)— EATC3 , 1  )>*VCN)  +  EATC3 , 2)*PCN)+EATC3 , 3)*PHICN)+EATC3 , 1) 

A-RCN) 

2100  CONTINUE 

. . . . .  TIME  RESPONSE  i  «"■  » ■  ■  *" . . . .  n '»■»’  ■  **■> 

. **°**° . °  nwwxt N+l ) -XCN+1 )  ♦  CONTROL  VALUE  »"■’ *»■ . mmi  imn*— * 

2200  VCM)«VCM)+C0NC1) 

2210  PCM)»PCM)*C0NC2) 

2220  PHICM)-PHICM)+C0NC3) 

2222  ARC— VC  M)  /VEL 

2225  BETACM)-ASINCARG) 

2230  RCM)»RCM)+C0NC4) 

f  Cm  ■■■  i"»  **°**  I'm i INCREMENTAL  TIME  RESPONSE  OUTPUT**  tmmtmummmmmm  mm 
<*  9050  VRITEC1 , 9060)  T,  VCM)  ,  PCM)  ,  PHICM)  ,  RCM)  ,  BETACM) 
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9060  F0RMAT<2X,F5.2,1X,F8.3,3X,F7.3,3X,F7.4,4X,F6. 3,4X,F7.5> 

9070  N*N+1 

9080  IF  CN.LT.ITD  GOTO  1510 

9090  CONTINUE 

9100  END 

SINCLUDE: 'GMINV. FOR' 

SINCLUDE: ' VADD. FOR' 
aiNCLUDE: 'DOT. FOR' 
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F-15A  AIRCRAFT  DESCRIPTION 


The  following  data  and  diagram  were  extracted  from  AFFTC 
TR— 7 6-48  CIS]. 


Horizontal  Stabilators: 
Area  (Total)  _____ 

MAC _ 

Span _ 


111.36  sq  ft 
3.28.  ft 
28.25  ft 


Vertical  Tails  (Each) : 


Area 

MAC 


Spaa  (Exposed) 


62.61  sq  ft 
6.75  ft 
10.32  ft 


OCtmPL  SURFACE  AREAS 


Flaps  (Tbtal)  _ 

Ailerons  (Total) 
Rudders  (Total) 
Speedbralce 


35.84  sq  ft 
26.48  sq  ft 
19.94  sq  ft 
31.50  sq  ft 


Figure  El.  F-1S  General  Aircraft  Layout 


121 


F-15A  Aircraft  Data 


Aircraft 

Length 

Height 

Wetted  Area 

Takeoff  Gross  Weight 

Internal  Fuel  Capacity 

Wing 

Reference  Area 
Span 

Aspect  Ratio 
Taper  Ratio 
Incidence 

Leading  Edge  Sweep 
Mean  Aerodynamic  Chord 
Theoretical  Root  Chord 

Horizontal  Stabilators 

Area 

Span 

Deflection  Limits 
Mean  Aerodynamic  Chord 
Root  Chord 

Tc 1  Length  (.25  c  to  .  25 

V 


62.52  Ft 
18.63  Ft 
2668. 8  Ft 
39,770  Lbs 
11,138  Lbs 


608  Ft2 
42.81  Ft 
3.0 
0.25 
0  Deg 
45  Deg 
15.94  Ft 
22.77  Ft 

120  Ft 
28.25  Ft 
♦15  to  -29  Deg 
8.27  Ft 
11.43  Ft 
T)  20.08  Ft 
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Vertical  Stabilizers 


Area 

Span 

Mean  Aerodynamic  Chord 
Root  Chord 

Tail  Length  C .  25  c  to  .25  c 

V  V 


125.22  Ft 
10. 32  Ft 
6.75  Ft 
9.58  Ft 


17.69  Ft 


KEY  INS/AIR  DATA  OUTPUTS 


INERTIAL  MEASURING  UNIT 
CN-137G 'ASN-109 


;  j  I 

i  ii 

\'ii  ii 


a  c 


PHYSICAL  CHARACTERISTICS 

-  Weight  -  36.5  pounds. 

-  Volume  -  1658  cubic  inches. 


IMU  EQUIPMENT  RACK 
MT-4390/ ASN-109 


-  Power  Required  -  115v,  3  phase,  400  Hz,  2S5  watts, 


Figure  FI.  F-15A  HIDEC  INS 


A!3  2AYA  COMPUTES 
(SPERRY) 


AIR  DATA  COMPUTER 
AN/ASK-6 


CHARACTERISTICS 

-  Digital  Computer. 

-  Continuous  BIT  and  Initiated  BIT  capability  provided. 

-  Initiated  BIT  available  on  ground  ONLY. 

-  Weight  -  12  pounds. 

-  Volume  -  508  cubic  inches. 

-  Power  Required  -  11 5v,  400  Hz,  53  watts. 


Figure  F2. 


F—I5A  HIDEC  Air  Data  Computer  System 


TABLE  FI 

tCEY  SENSOR  OUTPUTS 


Quantity 

Source 

0 

INS 

<P 

INS 

INS 

p 

Flight  Control 
Computer 

q 

Flight  Control 
Computer 

r 

Flight  Control 
Computer 

a 

X 

INS 

a 

y 

INS 

a_ 

Flight  Control 
Computer 

VE 

INS 

VN 

INS 

VDVN 

INS 

VT 

ADC 

h 

INS/ADC 

Mach 

ADC 

P 

ADC 

q 

Computed 

w 

Computed  from  Fuel 
Totalizer  and  Basic 
Weight 
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F-15  STABILITY  DERIVATIVE  MODELS 

The  following  figures  depict,  the  developmental  flight 
test  data  used  to  formulate  stability  derivative  models  for 
the  F-15A  Inflight  a  Estimator.  In  all  cases,  flight  test 
data  precedes  modeled  parameters  used  in  the  estimator. 
References  for  documents  from  which  the  actual  flight  test 
data  plots  were  extracted  from  are  included  immediately 
following  the  figure  title.  All  data  and  graphs  are 
unclassified. 


Mo^iht  or  Itucw 
(W.O«* FT*  *•<>'*) 


GROSS  WEIGHT  (LBS) 

Figure  G6.  F-iSA  Estimator  IXX  Model 
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CvO  l  X  Cold  -  S9I11S  )  A  A I 


I 


» 


F-15A  ANGLE  OF  ATTACK  ESTIMATOR 


IYY  PROGRAM  CODE  MODEL 
SECOND  ORDER  POLYNOMIAL  FIT 


GROSS  WEIGHT  (LBS) 


Figure  G7.  F-iSA  Estimator  IYY  Model 
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IZZ  (SLUGS  -  FT2)  X  I0A3 


APPENDIX  H 
ELEXI  PROGRAM  AO A 
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ELEXI  PROGRAM  AO A 


The  program  AO A  is  the  modification  of  the  basic  inflight 
a  estimator,  altered  for  implementation  on  the  NASA  ELEXI 
computer.  The  computer  allows  both  past  flight  and  actual 
real  time  inflight  data  reduction.  The  enumerated  input 
signals  are  read  in  and  the  program  then  sequentially  reduces 
each  time  segment  of  data  before  moving  on  to  the  next 
increment.  The  result  is  similar  to  onboard  implementation  of 
the  algorithm  for  real  time  estimation  of  angle  of  attack. 


PROGRAM  AOA 


THIS  PROGRAM  UTILIZES  AIRCRAFT  GEOMETRY.  6  STABILITY  DERIVATIVES.  AND 

inputs  from  the  ins  ano  caoc  to  calculate  angle  of  ATTACK  OF  THE 

AIRCRAFT  IN  ANY  FLIGHT  CONDITION.  (VERS'ON  2  CLEAN  CONFIGURATION). 


INPUTS  FROM  FLIGHT  OATA  FORMAT: 

FUEL  -  FUEL  WT  ON  BOARO  (LBS) 

-  currently  a  calculated  parameter 

THETA  -  PITCH  ANGLE 

-  enC3 ,  ADC/ I  NS  (DEG) 

PHI  -  BANK  ANGLE 

-  eno4,  ADC/ INS  (OEG) 

PS  I  -  HEADING  ANGLE 

-  eno7 ,  ADC/ I  NS  (DEG) 

P  -  ROLL  RATE 

-  dd25,  OFCC  (DEG/SEC) 

Q  -  PITCH  RATE 

-  ddOl .  OFCC  (DEG/SEC) 

fl  -  YAW  RATE 

-  <3028 ,  OFCC  (OEG/SEC) 

AX  -  BODY  X  ACCELERATION 

-  IbSD.  INDICATED  LONGITUDINAL  ACCELERATION  AT  C.G.  (G) 

AY  -  BODY  Y  ACCELERATION 

-  dd27 ,  LATERAL  ACCEL  <FT/SEC*»2) 

AZ  -  BOOY  Z  ACCELERATION 

-  dd02.  NORMAL  ACCEL  OFCC  '(DEG/SEC) 

VE  -  VELOCITY  EAST 

-  enOS,  INS  (FT/SEC) 

VN  -  VELOCITY  NORTH 

-  enOS ,  INS  (FT/SEC) 

VOWN  -  VELOCITY  OOWN 

-  enlO,  INS  (FT/SEC) 

VT  -  TRUE  AIRSPEED 

-  I g22 ,  AOC/INS  ( KTS ) 

H  -  ALT  I TUOE 

-  en09,  ADC/ INS  (FT) 

RMACH  -  MACH  NUMBER 

-  I g23 .  ADC/ INS 

RHO  -  AIR  DENSITY  ( ) 

-  1025  •  s.l.a.  density 

COOT  -  PITCH  ACCELERATION.  OFF  INS  (DEG/SEC)  ••  ESTIMATED 

-  ESTIMATED  FROM  Q 


EXTERNAL  OPENR ,  OPENW.  FREAO ,  FWRITE,  CLOSER.  CLOSEW 
LOGICAL  OPENR.  OPENW.  FREAO,  FWRITE 
CHARACTER* 1 8  3IGS(50) 

REAL'S  ROATAOO).  N 

INTEGER'4  UNIT).  UNI  TO,  NIN.  NOUT.  NAVA  I L 


PARAMETER  (  UNIT  I  -  3.  UN  I  TO  -  7,  NIN  -  18,  NOUT  -  29) 

OATA  SIGS  /  'FUEL',  ’THETA’,  'PHI’,  ' PS  I ' .  ’P',  ’O'. 

'R'.  ’AX',  ’AY*.  ’AZ'.  ’VE’,  ’VN'. 

‘VOWN’,  ’VT’,  '  H  ’  ,  ’  RMACH '  ,  ’RHO’,  ’  A  I  NF  ’  , 

‘  AWB  ’  .  ’VX’,  ’VV,  ‘OIR’ ,  ’ALPHAG'  . 

’  VWN  ’  ,  '  VWE  '  ,  ’VWZ‘,  'QOOT  ‘  ,  ’VX’.  ‘VZ’, 

21  *  '  '  / 


•  OPEN  OATA  FILES  AND  SET  UP  OUTPUT  FILES 

« 

I F (  .NOT.  OPENR ( UN  I T I  ,  ’ FLDT .IN', NAVA  I L ) )  THEN 
Wfl I TE ( 7 , 100) 

STOP 
ENO  IF 

* 

IF<  .NOT.  OPENW(UN  I  TO , ’ FLDT .OUT ’ , NOUT . S I GS , ‘CMP2 ' ) )  THEN 
WRITE (7, 200) 

STOP 
ENO  IF 

100  FORMAT ( // 'OOPS  1’/) 

200  FORMAT(//’OOPS  2'/) 


*  INITIALIZE 

ft 

PROGRAM  PARAMETERS 

ft 

OL 

Mi 

0.0 

1000 

S 

m 

603.0 

1  100 

XT 

m 

17.7 

1  1  10 

XSM 

■» 

0.0300  •  13.94 

1200 

CBAR 

m 

13.94 

1300 

S 

m 

42.81 

MOO 

WTORY 

m 

29303.00 

1300 

ft 

STORES 

m 

0.0 

PI 

m 

3. 141392834 

OTOR 

m 

PI  /  180.0 

•  REAO  IN  FLIGHT  OATA  PARAMETERS  FOR  SPEC iF 1C  TIME  INCREMENT. 
300  I F (  FREAO(UNIT)  , T I  ME , ROATA ) )  THEN 


FUEL 

-  ROATA( 

1  ) 

THETA 

-  ROATA ( 

2) 

PHI 

-  ROATA< 

3) 

PS  1 

•  ROATA< 

4) 

P 

-  ROATA ( 

3) 

Q 

-  ROATA ( 

0) 

R 

•  ROATA( 

7) 

AX 

-  ROATA* 

8)  •  32.1/4 

AY 

-  RDAT. .  ( 

3) 

AZ 

-  -1.0  • 

RDATA( 10) 

VE 

-  ROATA ( 

"1) 

VN 

-  ROATA( 12) 

VDWN 

•  RDAT A ( 13) 

VT 

-  ROATA( 14)  *  1 .687778 

H 

-  ROATAC 13) 

RMACH 

•  ROATA(  M) 

RHO 

•  ROATA » 

)  •  0.0023769 

A  INF 

•  ROATA( 18) 

WT 

•  WTORY 

♦  STORES  ♦  FUEL 

K  • 


CHANGE  ANGULAR  PARAMETERS  FROM  DEGREES  TO  RADIANS  ANO  CALCULATE 
TRIGOMETRIC  FUNCTIONS.  CALCULATE  QOOT. 

PHI  -  PH!  •  OTOR  r 

PS  I  •  PS  I  ■  CTOR 

THETA  -  THETA  •  OTOR 

P  -  P  •  OTOP 

Q  -  Q  •  OTOR 

R  -  R  •  OTOR 


COSPHI  -  COS ( PH  I  ) 
COSPSI  -  COS ( PS ! ) 
COSTHA  .  COS< THETA) 
SINPHI  -  SIN(PHI ) 
SINPSI  -  S I N( PS  I ) 
SINTHA  -  SIN(THETA) 

COOT  -  (Q  -  QL)  *  20.0 
0L  -  Q 


•  MOOEL  STA8ILITY  DERIVATIVES 
WTS  -  WT / 1 ,OE*OS 

4200  R i XX  -  1 . 1 1 9E *04  •  (WTS**2)  ♦  4.087E*03  •  WTS  ♦  2.025E*04 

RIYY  «  7 . 428E  *03  •  (WTS**2)  ♦  8.723E*03  •  WTS  ♦  1.502E*0S 
R 1 22  -  1.454E*04  •  (WTS**2 )  ♦  S.908E+03  *  WTS  *  1.878E+03 
CMO  -  0.018 


*  ESTIMATE  VX  ANO  VZ 


VX  -  (  COSTHA  •  COSPSI  )  •  VN  ♦ 

(  COSTHA  •  SINPSI  )  *  VE  - 
(  SINTHA  )  •  VDWN 
VX  «  ABS(VX) 

VZ  -  (  (SINPHI  •  SiNPSI)  ♦  (COSPHI  •  SINTHA  •  COSPSI)  )  •  VN  < 
(  (-SINPHI  •  COSPSI)  ♦  (COSPHI  •  SINTHA  ■  SINPSI)  )  •  VE 
(  COSPHI  •  COSTHA  )  •  VOWN 


•  COMPUTE  NO  WIND  ESTIMATE  OF  ALPHA  BASED  ON  INERTIAL  VELOCITIES  IN  THE 

•  BODY  DIRECTIONS. 


VZABS  -  ABS(VZ) 

ALPHAG  -  ATAN( VZA8S/VX ) 

IF  (VZ  .LT.  0.0)  ALPHAG  -  -ALPHAG 


•  DERIVE  LOAO  FACTOR  FROM  AZ 


8000  N  -  (  -1.0  •  (  COS ( ALPHAG )  •  AZ  ♦  S I N( ALPHAG )  •  AX  )  l  32.2  ) 

♦  1.0 


•  CALCULATION  OF  ALPHA  ESTIMATED 


CL  -  (  (N*WT*XT)  ♦  (QDOT«RIYY)  ♦  ( P*R* ( R I XX-RI ZZ )  )  * 

CMO*  ( 0 . 5*RH0 * ( VT**2 )  )  *S*CBl.R  )  / 

(  (0 . 5“RHO* ( VT**2 ) )  •  S  *  XT  •  (  1 . 0* ( XSM/XT ) )  ) 
AWB  -  0.79  ♦  1S.44*CL  -  2.75«RMACH 
AWBR  -  AWB  •  OTOR 


D  WIND  CALCULATION 

COS AWB  -  COS (AWBR) 

SINAWB  -  SIN(AWBR) 

VAN  -  (  COSPSI  ■  COSTHA  •  COSAWB  ♦ 

COSPSI  •  SINTHA  ■  COSPHI  •  SINAWB  ♦ 

SINPSI  •  SINPHI  •  SINAWB  )  •  VT 

VAE  -  (  SINPSI  •  COSTHA  •  COSAWB  ♦ 

SINPSI  *  SINTHA  •  COSPHI  •  SINAWB  - 

COSPSI  •  SINPHI  •  SINAWB  )  •  VT 


VAZ  -  (  -1.0  •  SINTHA  •  COSAWB  ♦ 

COSTHA  •  COSPHI  •  SINAWB  )  »  VT 


1X2 


VWN  •  VN 


VAN 


VWE  m  VE  -  VAE 


VW2  -  VOWN  -  VA2 


•  STANDARO  WIND  CALCULATIONS 

VW  -  SQRT(VWN«*2  +  VWE«*2) 

ARG  -  VWN/VT 

0  I R  -  ACOS  <  ARQ )  ♦  3.1714 

IF  (VWE  .LT.  O.O)  DIR  -  OIR  -  3.174 

VV  -  -1 . 0*VW2 


•  SET  UP  NEW  OUTPUT  PARAMETERS  ANO  WRITE  OUTPUT  FILE. 

* 

ROATA( 19)  -  AW0 

RDATA ( 20 )  -  VW 

ROATA ( 2 1 )  -  VV 

RDATA( 22 )  -  DIR 

ROATA (23)  -  ALPHAG  /  OTOR 

RDATA( 24 )  -  VWN 

RDATA(2S)  -  VWE 

RDATA( 28 )  -  VWZ 

ROATA (27)  -  QOOr 

ROATA ( 28 )  -  VX 

R0ATA( 29 )  -  V2 

• 

CALL  FWR I TE ( UN  I  TO . TIME . ROATA ) 

« 

GO  TO  300 
ENO  IF 


CALL  CLOSER (UN I T I ) 
CALL  CLOSEW ( UN  I  TO ) 


STOP 

ENO 
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APPENDIX  I 


FLIGHT  TEST  CARD 


FLIGHT  TEST  CARD 


The  Following  is  the  actual  NASA  flight  test  card  used  for  the 
robustness  test  on  the  NASA  F-15A  HIDEC  aircraft.  Test  point 
22  is  the  modified  split-S  maneuver  designed  specifically  to 
test  the  a  estimator.  The  test  point  was  flown  once  as 
printed  and  repeated  once,  substituting  a  30  degree  left  wing 
down  bank  at  2  g’s  for  step  22B. 


F-15A  NO.  8 
FIT  NO.  515 


NASA  335 


16. 

17. 

18. 

19. 

20. 

21. 

22. 

23. 


Figure  II 


PAGE  6  of  6 
DATE 


20K70.5  -  0.9M 


MIL  THRUST  ACCEL.  0.5  -  0.9M  COUPLED 
(  228  -  423  KCAS  ) 

MIL  THRUST  5  -  6a  TURN,  0.9M  COUPLED  - 
10  SECONDS 

20K70.4  -  0.9M 


MAX  THRUST  ACCEL.  0.4  -  0.9M  UNCOUPLED 
(  182  -  423  KCAS  ) 

REPEAT  COUPLED 

20K70.6M  (  275  KCAS  ) 

MIL  THRUST  3  -  4g  TURN,  0.6M  UNCOUPLED  - 
30  SECONDS 

REPEAT  COUPLED 

20K70.7M  (  324  KCAS  ) 

PERFORM  THE  FOLLOWING  MODIFIED  SPLIT-S 
MANEUVER. 

A.  HDG  TRUE  NORTH 

B.  30  DEG  RWD  BANK  /  2g's 

C.  ROLL  INVERTED 

D.  MOD.  SPLIT-S  MANEUVER  /  4-5g's 

RTB 


NASA  F-15A  HIDEG  Flight  Test  Card 
Robustness  Test 
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